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ABSTRACT 



Dynamic positioning of a moving platform is analyzed using data from the 
Global Positioning System (GPS) acquired in Phase II of the Seafloor Benchmark 
Experiment on R'V Point Sur in August 1986. GPS position determinations are 
compared to simultaneous Mini- Ranger fixes. 

The GPS positions computed using only broadcast ephemeris data were within 
20 m from the Mini-Ranger fi.xes when data from four satellites were used and within 
30 m when data from three satellites and a geoidal height constraint were used. It was 
found that the position accuracy is degraded when data from a satellite reaching 
culmination is used. 
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THESIS DISCLAIMER 



The reader is cautioned that computer programs developed in this research may 
not have been tested for all cases of interest. While every effort has been made, within 
the time available, to ensure that the programs are free of computational and logical 
errors, they cannot be considered validated. Any application of these programs 
without additional verification is at the risk of the user. 
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I. INTRODUCTION 



The Seafloor Benchmark Experiment, a project of the Hydrographic Sciences 
Group of the Oceanography Department at the Naval Postgraduate School (NFS), was 
initiated in 1984 with a goal to establish stations on the seafloor in real-time with 
geodetic accuracy (Saxena, 1982). The first experiment (Phase I), conducted in .May 

1985, included a configuration of four benchmarks on the seafloor (Spielvogel et al., 
1987). Further study indicated a sizeable improvement in the positional accuracy of the 
inner stations as compared to the accuracy of the outer stations of a larger net 
configuration (Kumar and Saxena, 1985). For this reason, during the Phase II 
experiment in August 1986, a nine-station net configuration was used (Kumar, 1986). 
Three types of GPS receivers and a Mini- Ranger (MR) Falcon system were used to 
determine the position of the R/V Point Sur. 

Keeping the goal of establishing and continuing GPS research at NTS, Brown 
(1986) established the Texas Instruments 4100 GEOSTAR GPS software, written by 
the Naval Surface Weapons Center for their CDC Cyber 865 computer, on the NTS 
IB.M-3033 computer. Brown was able to validate the modified software using static 
position data, but not the KAL.MN2 program which uses dynamic position data. 

The main objective of this thesis is to determine the positions of a dynamic 
platform using GPS and to compare them with MR Falcon-determined positions. 
Other objectives were to install software for processing the raw data from cassettes to 
FIC (Floating-Integer-Character) format, modify the KALMN2 program and validate 
it using Phase II data. Programs to plot and analyze the continuous ship positions by 
GPS and by MR Falcon were also written. 

The data used in this thesis w^ere acquired during Phase II of the Experiment. 
The GPS data are Umited to those collected by the TI-4100 receiver from N'SWC 
installed on board of the R/V Point Sur. The day chosen for analysis w^as 16 August 

1986, since adequate pitch and roll data are available for that day. 
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II. INSTRUMENTATION 



A. POSITIONING SYSTEMS 

1. GPS receivers 

During the experiment, the following GPS receivers belonging to several 
agencies and companies were installed on the R/V Point Sur: 

- TI-4100 receiver from Naval Surface Weapons Centers (NSWC) 

- TI-4100 receiver from National Ocean Survey (NOS) 

- Eagle Mini-Ranger receiver from Motorola, Inc. 

- Trimble 4000A from Trimble Navigation 

There were shore based GPS receivers at the following stations; 

Beach Lab 

- TI-4100 from Pacific Missile Test Center (PMTC) 

- Trimble 4000A from Trimble Navigation 
Perrier 

- Eagle Mini-Ranger from Motorola Inc. 

Dome 

- Magnavox Manpack from GPS Joint Program Office. 

2. Other positioning System 

For navigation and positioning on board were a MR Falcon system from 
Motorola, Inc., and a Loran-C receiver. 

B. OTHER SYSTEMS 

.An acoustic positioning system from Oceano Instruments to acquire acoustic 
ranges (Kuo, 1985) was installed on the R/V Point Sur; it also gives the position of the 
ship in relation to the bottom benchmarks (acoustic transponders) and can be 
integrated with a satellite receiver for computation of their positions using the GPS 
coordinates of the ship's acoustic transducer. 

In order to allow for corrections due to the pitch and roll of the ship, there w'as a 
data acquisition system, based on an HP9826 computer, to collect pitch, roll and 
heading data every second. Due to interfacing problems there was no possibility of 
logging heading data. The ship's Data Acquisition System (SDAS) collected, among 
other data, the ship's heading with a sampling period of about 19 s. 



11 



ro 



C ANTENNAS AND TRANSDUCER POSITIONS 

The antenna fef the MR Falcon was mounted on the mast of the ship in order to 
increase its range. The antenna of the TI-4100 from NOS was installed close to the 
.MR Falcon antenna. All the other GPS antennas were mounted on an elevated 
wooden table above the ship's laboratories and close to the transducer position (Figure 
. 1 ). 




Figure 2.1 Relative positions of antennas and transducer. 
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III. DATA PROCESSING 



A. COORDINATE SYSTEM 
I. Parameters 

In order to compare the positional accuracies of the TI-4100 GPS receiver 
with the Mini-Ranger (MR) Falcon, computations were made in a single coordinate 
system. 

Although there was a Doppler survey to find the WGS 72 coordinates of all 
control stations used, due to the impossibility of setting a station over Dome Ecc, the 
site for one of the MR Falcon shore stations, the WGS 72 coordinates of nearby 
station Dome, 68 m away, were obtained. During the Phase II experiment, GPS 
provided coordinates in the WGS 72 system (Bomford, 1980, Sepplin, 1974), although 
it converted to WGS 84 on 1 January 1987 (Decker, 1986). 

The coordinates of Dome were computed in NAD 27 and then transferred to 
Dome Ecc. The coordinates of Dome Ecc were converted from NAD 27 to WGS 72 in 
order to compute directly the ship positions in this system, avoiding the problem of 
converting each position individually. 

The subroutines used for converting geodetic coordinates from WGS 72 to 
Universal Transverse Mercator (UTM) and vice-versa wer written in Fortran and are 
based on Basic subroutines used by DMA for an HP9826 desktop computer. A cross 
check was made which revealed a maximum conversion error of 0.01 m. which is well 
within the accuracy limits of the system and, hence, does not affect the analysis of the 
data processed. 

The following parameters were used: 

WGS72 

Semi major axis ; 6,378,135 m 

Flattening : 1/298.26 

UT.M projection 

Central meridian: 123® 00' 00" W 

If this package of subroutines is used for any location other than the area of the 
Seafloor Benchmark Phase II, the UTM parameters must be changed in order to get 
the correct conversions. 
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GPS receivers: Ferrler, Dome Ecc 

Beach Lab, Ferrler, Dome 



Depths In meters 
Distances in kilometers 



Figure 3.1 Area of the experiment. 

2. Coordinates of stations 

The two shore stations were located on high points along the coast (Figure 
3.1). Their locations and coordinates in WGS 72 are as follows; 
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Ferrier 

Latitude ^36° 33' 53".748 N 
Longitude” 121° 53' 52". 939 W 
Dome Ecc 

Latitude 36° 18' 20".897 N 
Longitude 121° 54' 00".594 W 

The elevations of the different antennas were set in the Mini-Ranger receiver and the 
measured slant distances were converted 

to appropriate horizontal distances by the receiver system during data acquisition. 

B. REDUCTIONS OF POSITIONS TO THE TRANSDUCER 
1. Offsets of transducer relative to antennas 

Since the two antennas, namely the .MR Falcon antenna near the top of the 
mast and the TI-4100 antenna on the elevated table above the lab on the R, V Point 
Sur, were at different locations, a common point had to be chosen to compare both 
systems. The acoustic transducer mounted beneath the hull was chosen as a common 
point, which was used for further integration with the acoustic data, to compute the 
position of the Seafloor Benchmarks. 

For data reduction there was a need for pitch, roll and ship's heading data. 
While the existing pitch and roll data were sampled with a period of 1 s, the ship's 
heading was sampled with a period of about 19 s. Computed headings were tried to 
improve the rate but tests indicated that the real data, even with a low sampling rate 
was better than the computed. The method used to compute the heading is described in 
Section C.3. 

The geometric offsets between the MR Falcon, TI-4100 antenna and the 
acoustic underwater transducer were computed using a right-handed coordinate system 
centered in each antenna, and having the Y axis oriented towards the bow, the X axis 
to starboard and the Z axis towards the zenith. Using this convention, all offsets of 
the transducer are negative. Using the formulas developed in the next section, the 
corrections for pitch, roll and heading data were computed and applied. 

For the MR Falcon positions only horizontal corrections were applied, since 
the positions are computed in a two-dimensional system. For the GPS positions the 
corrections were applied to all three coordinates. 
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For both types of positions the coordinates are converted to UTM, corrected 
and converted badk to WGS 72. For the GPS positions the Z correction is applied to 
the geoidal height. 

The following offsets were measured and computed by Prof. Tucker and Mr 
James Cherry, and changed to the above format. 

Transducer offsets relative to the GPS antenna : 

XOFF = -4.268 m 
YOFF = -2.955 m 
ZOFF = -9.505 m 

Transducer offsets relative to the MR FALCON antenna : 

XOFF = -5.345 m 
YOFF = -9.453 m 
ZOFF = -16.152 m 

The course was computed in a range 0® to 360®, and the pitch and roll data 
were stored in a range -90® to -1-90®; following this convention the roll is positive when 
the port is up and the pitch positive when the bow is down. The sign of roll was 
changed to conform with the mathematical model described next. 

2. Math model 

The method to compute the correction to the coordinates is based on a seven- 
parameter transformation described in Moffitt and Mikhail (1980), so successive 
rotations of the vector defined by XOFF, YOFF and ZOFF, give the corrections to the 
coordinates illustrated in Figure 3.2. 

The rotation matrices are defined as follows: 




1 0 0 

0 cos(PITCH) sin(PITCH) 

0 -sin(PITCH) cos(PITCH) 



M 



r 




0 -sin(ROLL) 

1 0 

0 cos(ROLL) 



cos(ROLL) 

0 

sin(ROLL) 
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A -Antenna T - Transponder 



Figure 3.2 Rotation of axis. 



M 



h “ 



cos(HEADlNG) 

-sin(HEADING) 

0 



sin(HEADING) 

cos{HEADING) 

0 



0 

0 

1 



The orientation matrix is defined by the following equation: 



MQ=MpM^M|j (eqn3.1) 

and the corrections are computed by multiplying the orientation matrix by the offset 
vector, to get the following equations; 

DX = XOFF cos(ROLL) cos(ilEADlXG) + (eqn 3.2) 

YOFF {cos(PITCH) sin(HEADING) + 
sin(PITCH) sin(ROLL) cos(HEADlNG)) + 

ZOFF {sin{PITCH) sin(HEADING) - 
cos(PITCH) sin(ROLL) cos(lIEADING)} 
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DY = XOFF cos(ROLL) sin(HEADING) (-1.) + (eqn 3.3) 

YOFF {Cos(PITCH) cos(HEADING) - 
sin(PITCH) sin{ROLL) sin{HEADIXG)} + 

ZOFF {sin(PITCH) cos(HEADIXG) + 
cos(PITCH) sin{ROLL) sin(HEADIXG)} 

DZ = XOFF sin{ROLL) - (eqn 3.4) 

YOFF sin(PITCH) cos(ROLL) + 

ZOFF cos(PITCH) cos(ROLL) 

C. FALCON MINI-RANGER DATA PROCESSING 

1. Error of system 

The MR Falcon was used during the cruise for navigation and positioning. In 
post-processing it is used to evaluate GPS data. Since the Mini-Ranger is a system 
frequently tested and used, its error is known to be ± 3 m (Laurila, 1976; .Vlunson, 
1977). 

It would be better to have had more than two lines of position (LOP) to allow 
a minimization of the error in the computation of the Mini-Ranger position and also 
to allow a smoothing of the data in places where the data gets noisy as was seen in 
post processing. In Figure 3.3 it is possible to see the areas where the data is noisy and 
to see the small jumps characteristic of the system. 

2. Antenna positions 

For computing the antenna positions the coordinates of the shore stations 
were converted to WGS 72 UTM, and then each position was computed in the WGS 
72 UTM projection. The subroutine used to compute the positions was picked from 
program L'CO.MPS, a utility package used in the Hydrographic Sciences Group of the 
Oceanography Department. 

3. Heading Calculation 

As mentioned previously, the ship's heading was needed to compute the 
transducer's position, but there was no data available with the 1-s periodicity needed, 
so headings were calculated from Mini- Ranger data. The calculation of the ship's 
headings was made during the computation of the two antenna positions and kept on a 
separate hie. These computed headings were affected by the pitch and roll of the ship 
and also by the noisy data of the Mini-Ranger. 
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Figure 3.3 Plot of ship's track using MR Falcon data. 
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Some trials were made using different running averages for smoothing the 
computed headinf_in order to get a course close to reality. It was found that the 
courses obtained, even if filtered by a running average, oscillated and hence were not 
suitable for this analysis. 

4. Transducer positions 

The transducer positions were computed using the model described in 
Section B. The pitch and roll data used were sampled with a periodicity of 1 s most of 
the time when there was positioning data. The ship's heading was sampled with a 
periodicity of from 19 s to I minute. Although tests were made in order to compute the 
heading, as stated in the previous section it was decided to use recorded data, 
constraining the analysis to periods of constant course. 

5. Programs. Inputs and outputs 

During the development of the programs many more routines were written 
than mentioned here. In order to make the use of the programs easy, the total number 
of programs to process the .MR Falcon data was reduced to four. 

A unique format for the different position files was defined, allowing only one 
plotting program in which the user must change the title according to the type of 
positions plotted. 



Program: FALCON 

Input: .Month, day, year, hour, minutes, seconds (date time tag), 

code I, range I, strength I, codel, range2, strength2 
Output: Listing with date time tag and geographic positions 

File with date time tag and geographic positions 
File with date time tag and rough course 
Source: See Appendix A 



Program: COURSE SMOOTH 

Input; File with date time tag and rough course 

Output: File with date time tag and smooth course 

Listing with date time tag and smooth course 
Source; See Appendix B 
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Program: 

Input: 

Output: 



Source: 



TRANSDUC FALCON 

jFile date time tag and geographic positions of antenna 
Tile with date time tag and geographic positions of transducer 
Optional listing with date time tag and geographic positions 
of transducer 
See Appendix C 



Program: PLOT 

Input: File with date time tag and geographic positions 

Output: Plot with track of the ship in UTM projection 

Source: See Appendix D 



All files are in free format, and, if for by any reason a fixed format is needed, the user 
should change the programs to get the output desired. The MR Falcon data were 
originally on floppy disks in HP9816 format. A program was wxitten in Basic in order 
to convert and write the data onto a 9-track magnetic tape. The transfer of data to 
mass storage on the NPS IBM mainframe computer and data manipulation were made 
using existing procedures for the mainframe (Favorite, 1986; Mar, 1984). 

The scale and the area of the plot are fixed for the particular area of the 
Benchmark Experiment but may be changed in the program; thus any proportional size 
can be obtained by changing the control card as is e.xplained in the manual (W.R. 
Church Computer Center, 1981). 

.As the plot program can be used for plotting the outputs of different 
programs, the title should be changed accordingly, leaving the appropriate lines 
uncommented. An hour-minute tag is plotted every 5 minutes. 

D, GPS DATA PROCESSING 

1. Reading the TI-4100 cassettes 
a. System configuration 

The data from the TI-4100 receiver is recorded on digital cassettes in a 
format not directly compatible with the existing programs at NPS. To read and decode 
the data, it was necessary to purchase a MEMTEC cassette terminal (Model 5450XL) 
and obtain the appropriate softw’are from the Applied Research Laboratories (ARL) of 
the University of Texas at Austin. The MEMTEC cassette reader was connected to an 
IBM PC having a hard disk of 20 Mb (Megabytes), a modem and a math co-processor. 
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Preliminar>' tests made to verify the software and hardware connections of 
the tape deck witlLthe IBM PC indicated good results. During the tests it was found 
that the format of the final output of this package was not the same as defined 
previously by Scott and Peters (1983) but was in a new format, FIC, where each block 
is organized in a structured way with control records followed in order by all floating 
type data, integer type data and character type data. The first program, CON9TR, a 
part of the library installed by Brown (1986) in the SEF (Standard Exchange Format), 
was replaced by a new program using the FIC format. 

b. Reading cassettes to disk 

Each cassette is dumped to disk using the program MFERD. The operation 
takes about 20 minutes, and the file created fits on one floppy disk which can be kept 
for raw data backup. 

c. Converting the data 

The binary image of the cassette is then converted to a FIC format file 
using the program GS2FIC, which takes about 20 minutes. This format is suitable for 
processing in the IB.M-PC with appropriate software to be developed. 

In order to transfer these files to the main frame or other computer, it is 
necessary to run a program FICFICA that converts the FIC files to an ASCII, format 
where the data is organized in 80-column records. This operation takes about 30 
minutes, and the file created by this program can be dumped to a magnetic tape or 
sent to a host computer. 

The files created by these two programs cannot be stored on a normal 
floppy disk and were erased as soon as the data was on the mass storage of the 
mainframe and tested. 

d. Transfer to the mainframe 

The transfer of the files to the mainframe was done using a micro-computer 
connected via modem to the main frame. A terminal emulator distributed by the NFS 
Computer Center, SIMPC, simplifies the transfer of files (Simware Inc, 1984). 

Each file, corresponding to one cassette, was sent to a disk. This operation 
took about 3.5 hours, and sometimes the transfer was stopped in the middle. When this 
happened the transfer had to be started from the beginning. 

Each file has around 15000 lines and occupies about 65% of an A-disk. 
Using an extra disk every file was converted to a job format in order to transfer it to 
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mass storage. The data were sent to a member of a partition data set where all 
cassettes were stoi;ed. 

2. Conversion of formats 

After a study of the FIC format, the input of KALMN2 program, and the 
program CON9TR, a program was written called CVFICA, which had, after some 
debuging and improvements, a final form that fulfilled the needs for processing the 
data for this thesis. 

CVFIC.A, decodes each FICA block, producing two files and a listing of 
warnings and general information. One file contains the data to be processed by the 
KALMN2 program and the other has the positions computed by the GESAR software 
installed in the receiver. 

In the present version of the program the pseudo-ranges and other tracking 
information are stored only after the navigational data for the desired number of 
satellites is in the output file, since this data is needed to compute the positions in the 
Navigational Mode. It is very easy to change this program to process the data using 
the precise ephemeris. The program flags as bad the data that belong to a Space 
Vehicle (SV) whose data does not correspond to the existing SV identification in that 
tracer. The data are also flagged as bad when the status vector has any value other 
than zero. 

The program is able to convert the following data block types; 

101 - GESAR Versions 1.0+ Input data, 

3 - GESAR solutions, 

6 - Tracking data, 

8 - Tracking configuration, 

109 - Navigation VIessage record as transmited, 

9 - De-blocked subframes 1-3 from block 109 
1 1 - Receiver error block, 

13 - Tape header,' trailer. 

According to the CVFICA program any other block is read and dumped to 
the listing with a warning saying that the program was not able to handle it. No such 
warning was found after using the program with data from 10 cassettes. 

The hsting also indicates when the navigational data for each satelUte is 
received in order to give to the user the status of the constellation of satellites. 
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As this program outputs the GESAR position solutions and the data to the 
KALMN2 prografn into data set files, these must be created prior to running the 
program. 

A control file is used to pass to the program information such as two lines for 
the title of the listing and the starting date of the GPS week. The starting date allows 
the conversion of the time tag in seconds to a date time tag. The program does not 
work with data that crosses the one-week limit. 

3. Program information 

Program: CVFICA 

Input: FICA files without the first two comment lines 

Control file with the title and starting date of GPS week 
(Appendix J) 

Output: File with GESAR solutions 

File in XSWC format (input to KALMN program) 

Listing with warnings and general information 

Source: See Appendix E 

4. Computation of positions 
a. Problems found 

When starting to use the Kalman filter program (Brown, 19S6) for 
processing the GPS data some difficulties were encountered, either due to the 
replacement of the CON9TR program by an equivalent one or by errors existing in the 
program KALMN2. The first problems were caused by the different units used in the 
ARL and NSWC software which were a consequence of the lack of documentation 
during the development of the program CVFICA; these problems were easily found 
and solved. The last problems took a long time to find, not only due to the complexity 
of the program, but also due to the nonexistence of processed data to check the results. 
A big effort was made in order to correct all errors, but only when two cassettes were 
sent to NSWC and processed there were all errors removed from the processing. It was 
found that the last two major errors were due to causes external to the KALMN2 
program, one an error in the ARL software that affected the ambient temperature by a 
factor of 10 and caused an excessive tropospheric correction and the other a bad 
constant in the control file of the KAL.MN2 program. 
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b. Changes made 

In ad(iition to the corrections, some changes were made to improve the 
outputs. Now the program prints a date time tag just prior to all outputs labeled with 
time tag in seconds of GPS week. A title is printed as a header of the output listing. 
Other small changes were made to improve the code. .All changes were made carefully, 
and the programs w'ere throughly tested. 

c. Data processing 

After inserting data from ten TI-4100 cassettes into the mass storage, a 
control file for the KALMX2 program was made (Brown, 1986). For an initial 
approximate position, the position of the receiver based on the GESAR solutions with 
the time tag closest to the starting time was used for the data processing. For the 
receiver time bias, an iterative procedure was used; thus, the data were processed using 
a value, starting with zero, that was replaced by one listed in the outputs of the 
KAL.MN2 program when the computed positions were close to the predicted ones. 

Due to the good initial results of the test runs with certain satellite 
configurations, further tests were investigated for the different types of solutions. 

The value of -37 m was used as geoidal height in order to constrain the 
positions to the geoid when using data from less than four satellites or when forced by 
the control file. This value was obtained by interpolation on a chart (Blaha et al., 1986) 
and corrected for the antenna height above sea level, giving a corrected value of -31 m. 
The control file was set up in such a way as to save the positions on another file. The 
broadcast ephemeris was used, as the precise ephemeris w^as not available when the 
data were processed. 

d. Math model 

To compute position, the program uses the pseudo-ranges corrected for the 
various factors affecting them and an eight-state Kalman filter (Brown, 1986). 
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e. Program information 
Program : 'JCALMN2 



Input: 


File in NSWC format produced by the program CVFICA 
Control file as defined by Brown (1986)(Appendix I) 
Control file with the title and starting date of GPS week 


Output: 


Listing with positions, errors and other information 
as selected by the control file 
File with time tag in seconds of GPS week, X, Y, Z 
and Lat, Lon, HT, both in WGS 72 
and the number of satellites used in the solution 
Optional file with time tag 

and Geometric Dilution of Precision (GDOP) every minute 



5. Transducer positions. 

a. Procedure 

The positions of the transducer were computed using an algorithm similar 
to the one used with the MR Falcon data. The only difference is that the GPS 
positions are three-dimensional, so the corrections for all dimensions were computed 
and applied. 

The offsets between the transducer and the antenna and also the format of 
the input files were different. The program uses the same control file as the program 
CVFIC.'\ in order to compute the date time tag of each position. 

b. Program information 

Program: TRANSDUC GPS 



Input: 


File with time tag in seconds and WGS 72 coordinates 
Control file with the title and starting date of GPS week 


Output: 


File with date, time tag, latitude, longitude and ellipsoid height 
Optional listing with above information 


Source: 


See Appendix F 
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6. Plots 



It is possible to plot the track of the ship using either the GESAR solutions or 
the transducer positions as described before in Section C.5. Track plotting was done in 
order to give the GPS data coverage and data status. The program to plot is the same 
one used to plot the positions computed from the MR Falcon data but with the 
appropriate title (Figure 3.4). 
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Figure 3.4 Plot of ship's track using Tl-4100 data. 
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IV. DATA ANALYSIS 



A. COMPUTATION OF DIFFERENCES 

1. Approach 

To find the accuracy of the GPS positions in a two-dimensional coordinate 
system, the observed positions in WGS 72 coordinates are transformed to UTM 
coordinates and compared with the UTM coordinates of the corresponding MR Falcon 
positions interpolated for the same time tag. The differences in the coordinates of 
Mini- Ranger and GPS positions are expressed in meters. 

Due to the low sampling rate of the existing course data, only data from 
steady, straight courses were used for analysis. The data were divided into eleven 
periods of constant headings, while the data collected during the turns of the ship were 
deleted. This procedure made it possible to analyze the differences between GPS and 
.Mini- Ranger data during the different courses. 

Due to unexpected small differences found for positions computed using data 
from three satellites, several tests were made using various combinations of solutions. 
Although for data analysis many combinations were studied, three representative cases 
are discussed. 

Different starting times during processing of the GPS data using the 
KAL.MN2 program were used in order to see the effects of the propagation of errors 
due to the initial noisy data, but the discrepancies were inconsistent. The term "noisy 
data" is used to refer to the oscillatory behavior resulting from the Kalman filter when 
the number of satellites changes from three to four or vice versa. The effect of ship's 
motion on the transfer of antenna positions to two chosen common points, i.e. the 
acoustic transducer and the GPS antenna, was found to be 0.2 m, which falls within 
the accuracy limits of GPS and hence is considered negligible for this study. 

2. Case studies 

The following cases were analyzed; 

Case A : Data from four satellites were used or data from only three available 

satellites and geoidal height constraint for position computation were 
used. 
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Case B : Same as case A but using the geoidal height constraint when both 

three or four satellites were used. 

Case C : Same as B but not using data from SVll. 

3. Differences obtained 

a. Period I. 19:14 to 19:20 {Cases A and B) 



TABLE 1 

PERIOD 1, CASES A AND B. 

DIFFERENCES BETWEEN GPS AND MINI-RANGER POSITIONS 



Ranges! m) 

Avg Min- Max 

pos Svs ABA B 



179 


6,8,11 


26 


26 


24-32 


24-32 


80 


6.8,9,11 


24 


30 


21-31 


29-31 


1 


6,8.1^ 


27 


30 


27-27 


30-30 


1 


6,8.9.11 


34 


42 


34-34 


42-42 


1 


6.8.11 


34 


42 


34-34 


42-42 


23 


6,8,9.11 


27 


31 


25-29 


30-33 


1 


6,8,11 


28 


31 


28-28 


31-31 


74 


6,8,9,11 


23 


30 


20-30 


28-32 



* 

* 



* noisy data 



Case A - using four satellites, or 

four satellites with geoidal height constraint 
Case B - geoidal height constraint in case A 



pos - number of consecutive positions computed every second 



During the first part of this period, all positions were computed using data 
from three satellites. In both cases (Table 1) the differences found are in the same 
range, 24 to 32 m. For the next 80 positions, computed using data from four satellites. 
Case A lead to values lower than in Case B. The upper limit of the ranges in cases, 31 
m, is due to the initial oscillation from the Kalman filter. 

The next part of the period is one of noisy data, i.e. changes from four 
satellites to three and vice-versa. In both cases the positional accuracy is degraded, 
giving larger differences between GPS and .Vlini-Ranger positions. 

At the end of the period, with continuous data from four satelhtes and with 
the exception of one position, the differences are low, reaching 20 m in Case A and 28 
m in Case B. 
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It was found that the use of geoidal height as a constraint for position 
computations witb data from four satellites degrades the computed position. The 
reason for this is~ explained in Section 4.B. Noisy data causes oscillations in the 
observed differences which are due to the Kalman filter. 



TABLE 2 

SATELLITE VISIBILITY 
16 AUGUST 1986 



Time 


SV6 

El Az 


SV8 

El 


Az 


SV9 

El Az 


SVll 
El Az 


SV13 
El Az 


19:00 


59 


2 


51 


55 


24 


315 


57 


156 


— 


” 


19:20 


63 


21 


44 


44 


32 


321 


67 


147 


6 


198 


19:40 


66 


45 


36 


38 


39 


328 


75 


124 


15 


200 


20:00 


66 


72 


28 


34 


48 


335 


77 


78 


24 


203 


20:20 


62 


96 


19 


33 


57 


342 


71 


74 


33 


208 


20:40 


56 


115 


11 


34 


67 


350 


62 


36 


43 


214 



The offset of the GPS positions in relation to the MR Falcon positions is 
in a southeasterly direction (Figure 4.1 on page 40). Also the ship's course is close to 
the azimuth of SV9 (Table 2), which must be the reason for the successive tracking 
losses. The GPS antenna is in the shadow of the ship's mast when the satellite is low. 
b. Period 2. 19:20 to 19:25 {Cases A and B) 

All positions were computed using data from four satellites. There is a nine- 
second interruption of data in the beginning of the period corresponding to the time 
needed to change tapes in the receiver during logging. 

In Case A the differences range from 17 to 26 m (Table 3). The maximum 
value was found in the first position after the lack of data. Most of the differences are 
around 19 m. In Case B the differences range from 23 to 38 m. The minimum value is 
found after the gap in data. Most of the differences are around 29 m. Visual 
comparison is possible by referring to Figures 4.2 and 4.3. As before the offset of the 
GPS positions is in the southeasterly direction. 



31 



TABLE 3 

PERIOD 2, CASES A AND B. 

DIFFERENCES BETWEEN GPS AND MINT-RANGER POSITIONS 



pos SVs 



Ranges(m) 

Avg Min-Max 

ABA B 



295 6,8,9,11 19 29 17-26 23-38 



Case A - using four satellites, or 

four satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 

pos - number of consecutive positions computed every second 



c. Period 3. 19:25 to 19:30 ( Cases A and B ) 



TABLE 4 

PERIOD 3, CASES A AND B. 

DIFFERENCES BETWEEN GPS AND MINT-RANGER POSITIONS 







Av 


Ranges(m) 
g Min- Max 




pos 


Svs 


A 


'"B A 


B 




152 


6,8.9.11 


18 


30 17-19 


29-31 




1 


6 8: 11 


18 


32 18-18 


32-32 




1 


6.8,9,11 


13 


28 13-13 


28-28 




9 


6,8,11 


30 


31 17-37 


28-33 


’S'* 


118 


6 8 9,11 


26 


33 24-37 


32-36 




11 


6 8 11 


46 


36 31-49 


35-37 


❖ 


4 


6,8,9,11 


39 


37 36-43 


37-38 


❖ 


3 


6 8,ll 


37 


38 37-37 


38-38 





* noisy data 

Case A - using four satellites, or 

tnree satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 

pos - number of consecutive positions computed every second 



This period is characterized by having noise in different parts of the data 
and shows great variation in the differences. In Case A, at the beginning of the period 
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the average differences are IS m, after the first period of noisy data they reach values 
of 26 m, and at tf^e end of the period the differences are up to 46 m. In Case B at the 
begining of the peTiod the average differences are 30 m. This values changes to 33 m 
after the noisy data, and at end of the period the differences are up to 38 m. It is 
possible to distinguish the oscillation in the positions due to noisy data, much more 
emphasised in Case A than in Case B (Figures 4.4 and 4.5). The offset of the GPS 
positions is as before in the southeasterly direction. 

The successive periods lacking data from one satellite are as before due to 
SV9, which is still in an azimuth close to the ship's course. The fact that use of geoidal 
height as a constraint with data from four satellites for computation of positions does 
not improve the solution is evident in Table 4. 

d. Period 4. 19\30 to 19\35 (Cases A and B) 



TABLE 5 

PERIOD 4, CASES A AND B. 

DIFFERENCES BETWEEN GPS AND MINT-RANGER POSITIONS 







Av 


g 


Ranges^) 


■Ma.\ 




pos 


SVs 


A 


*B 


A' 


B 




7 


6 , 8.11 


38 


38 


38-40 


38-39 




40 


6 8 , 9.11 


31 


35 


29-40 


34-39 




7 


6 , 8.11 


41 


38 


32-43 


37-39 


❖ 


1 


6 . 8 , 9.11 


38 


37 


38-38 


37-37 




1 


6 , 8,11 


38 


38 


38-38 


38-38 


❖ 


49 


6 8 , 9,11 


30 


35 


27-38 


33-39 


❖ 


10 


6 , 8,11 


40 


37 


31-44 


35-38 


130 


6 , 8 , 9,11 


31 


38 


28-41 


37-40 




1 


6 , 8,11 


34 


38 


34-34 


38-38 




11 


6 , 8 : 9,11 


33 


39 


29-38 


36-45 


* 


9 


6 : 8,11 


41 


38 


33-46 


36-40 




2 


6 , 8 , 9,11 


39 


37 


38-40 


37-38 




1 


68 n 


38 


38 


38-38 


38-38 


55s 


81 


6 8 , 9,11 


38 


41 


35-44 


40-45 





noisy data 

Case A - using four satellites, or 

three satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A" 

pos - number of consecutive positions computed every second 
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This period is characterized by successive changes in the number of 
satellites used in the computation of positions; the differences are larger than observed 
in the previous periods.Successive losses of lock on SV9 were the cause for this. In 
Case A (Figure 4.6) it is possible to see the oscillations of the computed positions when 
the data are noisy. Improvement in the solutions is seen (Table 5) when we change 
from three to four satellites; the average differences go from 41 m to 30 m, becoming 
higher values when only three satellites are available. 

In Case B the differences have higher values than in Case A, but the ranges 
are smaller due to the smaller, noisy oscillations. As before the noisy is caused by the 
loss of lock on SV9, and the offset in the GPS positions is in the southeasterly 
direction from the Mini-Ranger positions. 

e. Period 5. 19:35 to 1940 (Cases A and B) 

This period has the noisiest data processed, so the differences have greater 
spread. Although the range between differences is wide, 20 to 68 m in Case A and 29 
to 67 m in Case B, the differences are mostly smaller in Case A than in Case B. 

In Figure 4.7 we see successive jumps from three satellites to four and vice 
versa. The GPS positions are south and east of the Mini-Ranger positions. The 
differences in the northings are almost double of the ones observed in the eastings. 
The reason for the frequent changes in the number of satellites used is the relatively 
low altitude of SV9 and an azimuth close to ship's heading, which together with the 
pitch and roll led frequently to shadowing of the TI-4100 antenna. 

/. Period 6. 19:47 to 19:53 (Cases A, B and C) 

This period is characterized also by noisy data, but it is for a different 
reason. Here data from SV8 is rejected during the phase of smoothing it \vith Doppler 
data. The observed roll is higher than in previous periods, caused perhaps by the 
almost 90° change in course. As Table 2 shows, SV8 is setting, having a lower 
elevation than any other satellite used during this period. 

There are no significant differences between Cases A and B, and average 
difierences range from 26 to 30 m. In Case C we see that even lowering the number of 
satellites to two or three and using the geoidal height as a constraint, smaller 
differences are obtained than the ones computed with data from SVll (Table 6). The 
offset in the positions has no relation to the course of the ship, i.e. it is still in the 
southeasterly direction, and a small or no time lag exists between the data of the two 
systems (Figure 4.8). 
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TABLE 6 



PERIOD 6, CASES A, B AND C. 

DIFFERENCES BETWEEN GPS AND MINI-RANGER POSITIONS 



pos 


SVs 


A 


Ranges(m) 

Avg 

B ®C A 


.Min- .Max 
B 


C 




16 


8,9.11.13 


28 


29 


29 


28-29 


28-30 


28-30 




1 


9,11,13 


27 


29 


28 


27-27 


29-29 


28-28 


* 


3 


8.9.flJ3 


28 


29 


28 


28-29 


29-30 


28-29 


»ic 


1 


9,11.1^ 


27 


28 


26 


27-27 


28-28 


26-26 




4 


8,9.11.13 


28 


28 


26 


27-28 


27-30 


25-27 


❖ 


1 


9.11,13 


28 


28 


26 


28-28 


28-28 


26-26 




79 


8,9.fl.l3 


27 


28 


26 


26-29 


26-30 


24-28 




1 




27 


28 


26 


27-27 


28-28 


26-26 


* 


35 


8,9,flJ3 


29 


30 


30 


27-32 


28-33 


27-33 




1 


9,11.1^ 


28 


29 


28 


28-28 


29-29 


28-28 


♦ 


17 


8,9.11,13 


29 


29 


29 


28-30 


28-30 


27-30 




1 


9,11,1^ 


28 


29 


27 


28-28 


29-29 


27-27 


♦ 


20 


8,9, f 1,1 3 


29 


29 


28 


28-31 


28-31 


26-32 




1 


9.11.1J 


28 


29 


27 


28-28 


29-29 


27-27 


♦ 


162 


8,9.11,13 


28 


29 


27 


26-31 


27-31 


24-31 




1 


9,11,1!5 


28 


30 


27 


28-28 


30-30 


27-27 




1 


8,9.fl,13 


29 


30 


28 


29-29 


30-30 


28-28 




1 


9,11,i!5 


27 


29 


27 


27-27 


29-29 


27-27 




1 


8,9.n.l3 


29 


30 


28 


29-29 


30-30 


28-28 




1 


9,1L13 


27 


28 


25 


27-27 


28-28 


25-25 




5 


8,9.fl.l3 


27 


27 


21 


26-27 


26-27 


20-23 




1 


9,11,13 


27 


28 


24 


27-27 


28-28 


24-24 


* 


2 


8,9.n.l3 


29 


30 


27 


29-30 


29-31 


20-21 




1 


9,11.13 


26 


27 


23 


26-26 


27-27 


23-23 




6 


8,9,11,13 


26 


26 


21 


26-28 


25-28 


20-24 





* noisy data 

Case A - using four satellites, or 

tnree satellites with geoidal height constraint 
Case B - aeoidal height constraint in Case A 
Case C - Ignoring SVl 1 but constraining 
with geoidal height 

pos - number of consecutive positions computed every second 



g. Period 7. 19:57 to 20:05 {Cases A, B and Q 

In Case B there is an improvement in the positions in relation to Case A. If 
we refer to previous periods, we see that what is happening is the reverse of what was 
happening there, i.e. the use of the geoidal height in positions computed using data 
from four satellites improves the solution. As before, we see that noisy data causes an 
oscillation in the observed differences but with small amplitudes. Noteworthy are the 
highly improved positions indicated by a 13-m difference (Table 7). 
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TABLE 7 



PERIOD 7, CASES A, B AND C 

DIFFERENCES BETWEEN GPS AND MINT-RANGER POSITIONS 



pos 


SVs 


A 


Ranges(m) 

Avg 

B C A 


Min-. Max 
B 


C 




185 


8,9,11.13 


25 


23 


21 


20-27 


19-25 


17-23 




1 


8 1L13 


24 


23 


20 


24-24 


23-23 


20-20 




22 


8,9,fl.l3 


22 


20 


18 


20-24 


18-22 


16-20 




1 


8,11,13 


22 


21 


18 


22-22 


21-21 


18-18 




17 


8,9,fl,13 


22 


21 


18 


20-25 


19-23 


17-21 




1 


8,11,1^ 


23 


22 


19 


23-23 


22-22 


19-19 




12 


8,9,flT3 


22 


20 


17 


21-22 


19-21 


17-21 




1 


8:11.1!} 


19 


19 


16 


19-19 


19-19 


16-16 


* 


72 


8,9,11.13 


20 


18 


16 


18-22 


17-19 


13-18 




1 


8,ll.l3 


19 


19 


16 


19-19 


19-19 


16-16 




13 


8,9,11.13 


19 


18 


14 


17-20 


17-19 


13-16 




1 


8,lf,13 


18 


17 


14 


18-18 


17-17 


14-14 


♦ 


81 


8,9.fl.l3 


21 


19 


16 


19-26 


17-24 


14-21 




1 


9,11.13 


17 


17 


14 


17-17 


17-17 


14-14 




35 


8.9.11.13 


20 


18 


16 


19-22 


16-20 


13-18 




1 


9.11,13 


17 


16 


13 


17-17 


16-16 


13-13 


* 


2 


8,9,fl,13 


29 


25 


23 


28-29 


25-26 


22-23 


* 


1 


9A\A^ 


27 


24 


21 


27-27 


24-24 


21-21 


* 


27 


8,9,11,13 


25 


23 


21 


18-36 


17-34 


15-31 





* noisy data 

Case A - usina four satellites, or . . 

three satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 
Case C - Ignoring SVII but constraining 
with geoidal height 

pos - number of consecutive positions computed every second 



The GPS positions are also south and east of the MR Falcon positions, but 
the fact of ignoring SVII causes a decrease in the differences in easting (Figure 4.8). 
h. Period 8. 20:07 to 20:13 {Cases A, B and C) 

During this period all differences in all cases are smaller than 21 m (Table 
8). The higher values are obtained after the computation of a position using fewer 
than four satellites in Cases A and B, and three satellites in Case C, and are due to the 
oscillation caused by the Kalman filter. Table 8 shows that constraining the positions 
computed using data from four satellites improves the solution; the minimum difference 
reaches only 10 m. Also, ignoring SVII causes a reduction of the differences, mainly in 
the easting, and we see that the differences are smaller than in the other cases (Table 
8). Figures 4.10 and 4.11 show the improvement in the positions when SVII is 
ignored. 
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=_ TABLE 8 

PERIOD 8. CASES A, B AND C 

DIFFERENCES BETWEEN GPS AND MINI-RANGER POSITIONS 



pos 


SVs 


A 


Ranges(m) 

Avg 

B C A 


.Min- Max 
B 


C 


37 


8,9.11J3 


19 


19 


14 


19-21 


16-18 


13-15 


1 


8.11. ll5 


18 


17 


14 


18-18 


17-17 


14-14 


153 


8,9,1 1J3 


20 


17 


15 


17-21 


16-19 


12-17 


1 


8:il.ll? 


17 


17 


13 


17-17 


17-17 


13-13 * 


28 


8.9.I1T3 


20 


17 


15 


19-21 


16-19 


13-16 


1 


8,ll,i;5 


17 


17 


13 


17-17 


17-17 


13-13 * 


142 


8,9,fl,13 


19 


16 


12 


17-20 


15-17 


10-14 



* noisy data 

Case A - using four satellites, or 

three satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 
Case C - Ignoring S V 1 1 but constraining 
with geoidal height 

pos - number of consecutive positions computed every second 



i. Period 9. 20:13 to 20:20 (Cases A, B and C) 



TABLE 9 

PERIOD 9, CASES A, B AND C. 

DIFFERENCES BETWEEN GPS AND MINI-RANGER POSITIONS 







Avg 


Ranges(m) 


Min- .Max 




pos SVs 


A 


B 


C 


A 


B 


C 


391 8,9,11,13 


20 


17 


13 


19-22 


15-18 


11-14 


1 no data 


18 


15 


11 


18-18 


15-15 


11-11 


31 8,9,11,13 

* noisy data 


20 


15 


11 


18-21 


15-16 


10-11 



Case A - using four satellites, or 

three satellites with geoidal height constramt 
Case B - geoidal height constraint in Case A 
Case C - Ignoring SVll but constraining 
with geoidal height . 

pos - number of consecutive positions computed every second 



37 



The differences found during this period are similar to the ones obtained in 
period 8. Table 9 mdicates that the differences in all cases are correlated as before. The 
oscillations are due to the noise in the data from the MR Falcon and not to the GPS 
positions, which define a smooth and continuous course (Figure 4.12). 
j. Period 10. 20:20 to 20:25 {Cases A, B and C) 



TABLE 10 

PERIOD 10, CASES A, B AND C. 

DIFFERENCES BETWEEN GPS AND MINT-RANGER POSITIONS 



pos 


SVs 


A 


Ranges! m) 

Avg 

B ^C A 


Min-Max 

B 


C 


9 


8,9,11,13 


20 


16 11 


19-21 


14-16 


9-11 


88 


8,9,13 


11 


11 11 


10-17 


10-15 


10-12 


206 


6,8,9,13 


25 


25 24 


15-26 


15-26 


15-26 



Case A - using four satellites, or 

tnree satellites \vith geoidal height constraint 
Case B - geoidal height constraint in Case A 
Case C - Ignoring SVll but constraining 
with geoidal height 

pos - number of consecutive positions computed every second 



During the first tests, the differences were smaller when the positions where 
computed from three satellites, SV 1 1 excluded and the ship constrained to the geoidal 
height, than the ones using four satellites, SVll or SV6 included. This happened 
because during data acquisition SV6 was selected to replace SVll, and while the 
receiver was not locked on the SV6, only data collected from three satellites was good 
(Figure 4.13). 

In Case B (Figure 4.14 or Table 10) constraint of the position to the 
geoidal height improves the solution where SVll is used. In Case C (Figure 4.15) we 
notice a reduction in the differences in the positions where data from SVll was used. 
When SV6 replaces SVll, in all cases the differences are larger but of same 
approximate magnitude. 
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k. Period II. 2(h30 to 20:37 {Cases A, B and Q 



TABLE 11 

PERIOD 11. CASES A, B AND C. 

DIFFERENCES BETWEEN GPS AND MINT-RANGER POSITIONS 



pos 

130 

287 



SVs 

6.8.9T3 

8,9,1^ 



Ranges! m) 

Avg ^ ' Min-Max 

A B ^C A B 



35 38 38 

36 36 36 



35-36 

30-37 



37-39 

30-38 



37-39 

30-38 



Case A - using four satellites, or 

tnree satellites with geoidal height constraint 
Case B - geoidal heiaht constraint in Case A 
Case C - Ignoring SVll but constraining 
with geoidal height 

pos - number of consecutive positions computed every second 



The initial positions were computed using data from four satellites, and 
after 20h32ml0.7s the receiver started to track SVI2 in place of SVS, but no broadcast 
ephemeris data were received during this period. The differences computed using data 
from three satellites are the same in all cases, but those computed using data from four 
satellites have different values. Thus in Case A the differences are smaller than in Cases 
B and C, where the positions are constrained to the geoidal height (Table II). In all 
cases the differences in eastings are larger than in northings, contrary to what was 
happening before; this offset is similar for the other two cases (Figure 4.16). 
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Figure 4.1 Period I (Case A). Using all available satellites. 
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Figure 4.2 Period 2 {Case A). Using all available satellites. 
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Figure 4.3 Period 2 (Case B). Constraining to geoidal height. 
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SERFLOOR BENCHMRRK EXPERIMENT 
PHASE 1 I 
16 AUGUST 1986 




Figure 4.4 Period 3 (Case A). Using all available satellites. 
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Figure 4.5 Period 3 (Case B). Constraining to geoidal height. 
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Figure 4.6 Period 4 (Case A). Using all available satellites. 
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Figure 4.7 Period 5 {Case A). Using all available satellites. 
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Figure 4.8 Period 6 (Case C). Ignoring SVll. 
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Figure 4. 9 Period 7 (Case B). Constraining to geoidal height. 
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SEAFLOOR BENCHMARK EXPERIMENT 
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16 AUGUST 1986 
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Figure 4.10 Period 8 (Case B). Constraining to geoidal height. 
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SEflFLGOR BENCHMnnK EXPERIMENT 
PURSE II 
16 RUGUST 1986 
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f MR FRLCON TRRNOUCER POSITIONS 

WDPs trrnoucer positions 
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Figure 4. 1 1 Period 8 (Case C). Ignoring SVl 1, 
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+ MR FflLCBN TflflNOUCEn POSITIONS 
KCP5 TRflNOUCER POSITIONS 
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Figure 4.12 Period 9 (Case C). Ignoring SVll. 
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SEAFLOOR BENCHMARK EXPERIMENT 
PHASE II 
16 AUGUST 1986 



2020 
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GPS TRRNOUCEfl POSITIONS 
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Figure 4. 1 3 Period 10 (Case A).. Using all available satellites. 
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JlfGPS TRRMOUCER POSITIONS 
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Figure 4.14 Period 10 (Case B). Constraining to geoidal height. 
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Figure 4.15 Period 10 (Case C). Ignoring SVll. 
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SERFLOOR BENCHMARK EXPERIMENT 
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16 AUGUST 1986 
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SV6, SV8, SV9, SV13 
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Figure 4.16 Period 1 1 (Case A). Using all available satellites. 



4. Programs. Inputs and outputs 

Two proems written for data analysis are summarized below; 



Program: COMPARE POSITION 

Input: Control file with period limits (Appendix K) 

File with MR Falcon positions 
File with GPS positions 

Output: Listing with date, time tag, and GPS geographic positions 

and differences in UTM coordinates 
Source; See Appendix G 



•Program: COMPARE PLOT 

Input; Control file with period limits (Appendix K) 

File with MR Falcon positions 
File with GPS positions 

Output; Plots of tracks using GPS and MR Falcon data 

Source; See Appendix H 

B. CONCLUSIONS 

The results of the data analysis using four satellites with the broadcast ephemeris 
indicate that the real-time positioning by GPS of a dynamic platform, e.g. a ship, under 
the best conditions is within 15 m but in most cases 20 m. 

Using data from three satellites and constraining the solution to the geoidal 
height, positional accuracies are within 20 m under the best conditions and 30 m in 
most of the cases. 

It was found that the GPS positions were always southeasterly of the Mini- 
Ranger positions. Due to lack of sufficient data and due to time limitations, correlation 
between the geometry of the satellites and ship's heading could not be investigated. 

The accuracy of the positions is correlated with the relative motion of the 
satellites; this was found when SVIl was ignored as in Case C, which improved the 
solution. If the variations in azimuth and elevation of all satellites are compared for all 
periods the data was processed, SVs 6 and 11 have variations larger than 110 degrees, 
while SVs 8, 9 and 13 have variations in azimuth smaller than 35 degrees (Table 2). It 
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can be seen that both SVs 6 and 1 1 reach their culmination during the period the data 
is analyzed. SV6 reaches its culmination around 20:15 and SVll around 19:50. 

When SVll^is approaching culmination, the solution is not improved when the 
position is constrained to the geoid; after culnaination however, application of the 
constraint improves the solution (Tables 12, 13 and 14). 

The degradation in the accuracy of the GPS positions using either SV6 or SVll 
is evident in Table 15 where the differences found in these situations are bigger than for 
positions computed without data from them. It is evident that even in Case C when 
SV6, with the same characteristics as SVll, is used to compute the position, the 
differences become larger. 

There are some theoretical studies (Landau, 1986) about the selection of satellite 
configuration in order to get the best results for computed positions. However, with 
the present data it is not possible to select a best combination, since the data 
acquisition w'as limited to a small observational period and a small number of satellites. 
Currently a way of measuring the effect of the geometry of the satellite configuration is 
through the analysis of the geometric dilution of precision (GDOP) parameters. These 
parameters include the position dilution of precision (PDOP), reflecting the dilution of 
precision in three dimensions, the the horizontal dilution of precision (HDOP), 
reflecting the dilution of precision in two dimensions, the the vertical dilution of 
precision (VDOP), reflecting the dilution of precision in the vertical dimension, and the 
time dilution of precision (TDOP), reflecting the dilution of precision in time (.Milliken 
and Zoller, 1980). Thus, a low PDOP provides a good geometric configuration. Table 
16 lists the observ’ed satellites and their GDOP parameters. 

It was found that the position of the antenna is important, and in the case of the 
NSWC TI-4100 receiver some error was caused by shadowing due to the ship's mast 
when the azimuth of the satellites w’as close to the course of the ship. 

The positions used for comparison with the GPS positions are themselves 
affected by the noisy signal of the MR Falcon positioning system. This problem could 
be reduced if more than two LOPs were used in the computation of positions. 
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TABLE 12 



MR FALCON POSITIONS VERSUS GPS POSITIONS. 
OBSERVED DIFFERENCES IN METERS DURING PERIOD 3 





Time 




Case 


A 




Case B 






Case C 


H 


M 


S 


DX 


DY 


DP 


DX 


DY 


DP 


DX 


DY 


DP 


19 


25 


0.7 


8 


-18 


19 


11 


-29 


31 


- 


- 


- 


19 


25 


5.7 


8 


-17 


19 


12 


-27 


30 


- 


- 


- 


19 


25 


10.7 


9 


-15 


18 


12 


-27 


29 


- 


- 


- 


19 


25 


15.7 


8 


-17 


19 


12 


-28 


30 


- 


- 


- 


19 


25 


20.7 


9 


-16 


18 


12 


-27 


30 


-- 


- 


- 


19 


25 


25.7 


9 


-16 


18 


12 


-27 


30 


- 


- 


- 


19 


25 


30.7 


9 


-15 


18 


12 


-27 


29 


- 


- 


- 


19 


25 


35.7 


9 


-16 


19 


12 


-27 


30 


- 


- 


- 


19 


25 


40.7 


10 


-15 


18 


13 


-26 


29 


- 


- 


- 


19 


25 


45.7 


10 


-16 


18 


13 


-26 


30 










DX, DY differences in the UTM coordinates 
DP derived position difference from DX and DY 

Case A - using four satellites, or 

three satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 
Case C - no data during this period 
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TABLE 13 














MR FALCON POSITIONS VERSUS GPS POSITIONS. 
OBSERVED DIFFERENCES IN METERS DURING PERIOD 6 


H 


Time 

M 


S 


DX 


Case A 
DY DP 


DX 


Case B 
DY 


DP 


DX 


Case C 
DY DP 


19 


52 


40.7 


17 


-22 


28 


19 


-22 


29 


16 


-21 


26 


19 


52 


41.7 


16 


-24 


28 


18 


-23 


29 


15 


-22 


27 


19 


52 


42.7 


16 


-24 


29 


18 


-23 


30 


16 


-22 


27 


19 


52 


43.7 


17 


-22 


28 


19 


-23 


30 


16 


-22 


27 


19 


52 


44.7 


15 


-24 


29 


16 


-25 


30 


14 


-24 


28 


19 


52 


45.7 


16 


-21 


27 


18 


-22 


29 


15 


-21 


26 


19 


52 


46.7 


10 


-27 


29 


12 


-28 


30 


8 


-27 


28 


19 


52 


47.7 


14 


-23 


27 


15 


-24 


28 


8 


-23 


25 


19 


52 


48.7 


14 


-23 


27 


15 


-23 


27 


9 


-20 


22 


19 


52 


49.7 


15 


-22 


26 


16 


-21 


26 


10 


-18 


21 


19 


52 


50.7 


17 


-20 


26 


18 


-19 


26 


12 


-16 


20 


19 


52 


51.7 


18 


-20 


27 


19 


-19 


27 


14 


-17 


21 


19 


52 


52.7 


19 


-19 


27 


20 


-19 


27 


15 


-17 


23 


19 


52 


53.7 


19 


-19 


27 


20 


-19 


28 


16 


-18 


24 


19 


52 


54.7 


13 


-27 


30 


13 


-27 


31 


10 


-26 


28 


19 


52 


55.7 


12 


-26 


29 


13 


-26 


29 


9 


-25 


26 


DX, DY differences in the UTM coordinates 
DP derived position difference from DX and DY 










Case A - using four satellites, or 

three satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 
Case C - Ignoring SVll but constraining with 
geoidal height 
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TABLE 14 

MR FALCON POSITIONS VERSUS GPS POSITIONS. 
OBSERVED DIFFERENCES IN METERS DURING PERIOD 8 



Case A Case B Case C 

DY DP DX DY DP DX DY DP 



Time 



H 


M 


S 


DX 


20 


7 


0.7 


11 


20 


7 


20.7 


10 


20 


7 


40.7 


10 


20 


8 


0.7 


9 


20 


8 


20.7 


10 


20 


8 


40.7 


10 


20 


9 


0.7 


11 


20 


9 


20.7 


11 


20 


9 


40.7 


11 


20 


10 


0.7 


ll 


20 


10 


20.7 


9 


20 


10 


40.7 


11 


20 


11 


0.7 


ll 


20 


11 


20.7 


12 


20 


11 


40.7 


12 


20 


12 


0.7 


12 



-15 


19 


9 


-14 


-17 


20 


8 


-15 


-15 


18 


9 


-14 


-18 


20 


7 


-17 


-17 


20 


8 


-16 


-17 


20 


8 


-16 


-17 


20 


8 


-15 


-17 


20 


8 


-15 


-17 


20 


8 


-15 


-15 


19 


9 


-14 


-17 


19 


8 


-16 


-17 


17 


9 


-14 


-16 


20 


8 


-15 


-16 


20 


9 


-14 


-16 


19 


8 


-14 


-15 


19 


9 


-13 



17 


5 


-13 


14 


17 


4 


-14 


14 


17 


4 


-12 


13 


18 


4 


-16 


16 


17 


5 


-14 


15 


17 


4 


-14 


15 


17 


5 


-14 


15 


17 


5 


-14 


14 


17 


5 


-14 


14 


16 


5 


-12 


13 


17 


5 


-14 


15 


16 


5 


-12 


12 


17 


4 


-13 


13 


17 


4 


-13 


13 


16 


4 


-12 


13 


16 


4 


-ll 


12 



DX, DY differences in the UTM coordinates 
DP derived position difference from DX and DY 

Case A - using four satellites, or 

three satellites with geoidal height constraint 
Case B - geoidal height constraint in Case A 
Case C - Ignoring SVll but constraining with 
geoidal height 
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TABLE 15 

MR FALCON POSITIONS VERSUS GPS POSITIONS. 
OBSERVED DIFFERENCES IN METERS DURING PERIOD 10 



H 


Time 

M 


S 


DX 


Case A 
DY DP 


DX 


Case B 
DY 


DP 


DX 


Case C 
DY 


DP 


20 


20 


7.7 


14 


-15 


20 


9 


-13 


15 


4 


-10 


11 


20 


20 


8.7 


14 


-15 


20 


9 


-13 


15 


4 


-10 


11 


20 


20 


9.7 


11 


-13 


17 


8 


-12 


15 


4 


-10 


11 


20 


20 


10.7 


9 


-13 


15 


7 


-13 


14 


3 


-11 


11 


20 


20 


11.7 


8 


-12 


14 


7 


-12 


14 


3 


-10 


11 


20 


20 


12.7 


7 


-12 


13 


6 


-12 


13 


3 


-10 


11 


20 


20 


44.7 


4 


-10 


10 


4 


-10 


10 


4 


-10 


10 


20 


20 


46.7 


4 


-10 


10 


4 


-10 


10 


4 


-10 


10 


20 


20 


47.7 




-11 


11 


3 


-11 


11 


3 


-11 


11 


20 


21 


35.7 


3 


-10 


11 


3 


-10 


11 


3 


-10 


10 


20 


21 


36.7 


3 


-10 


11 


3 


-10 


11 


3 


-10 


11 


20 


21 


37.7 


8 


-13 


15 


8 


-13 


15 


8 


-13 


15 


20 


21 


38.7 


10 


-14 


18 


10 


-14 


18 


10 


-14 


18 


20 


21 


39.7 


12 


-15 


19 


12 


-15 


20 


12 


-15 


20 


20 


21 


40.7 


13 


-16 


21 


13 


-16 


20 


13 


-16 


21 



DX, DY difierences in the UTM coordinates 
DP derived position difference from DX and DY 

Case A - using four satellites, or 

three satellites with geoidal height constraint 
Case B - geoidal heisht constraint in Case A 
Case C - Ignoring. SVll but constraining with 
geoidal height 
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TABLE 16 

PREDICTED GDOP PARAMETERS 
16 AUGUST 1986 



Time 




SV set 




GDOP 


PDOP 


HDOP 


VDOP 


TDOP 


19:00 


6 


8 


9 


11 


6.90 


5.75 


2.47 


5.19 


3.82 


19:20 


6 


8 


9 


11 


7.61 


6.30 


2.65 


5.72 


4.47 


19:40 


8 


9 11 


13 


3.37 


3.01 


1.97 


2.27 


1.53 


20:20 


8 


9 


11 


13 


3.86 


3.53 


2.75 


2.22 


1.55 




6 


8 


9 


13 


3.33 


2.94 


1.88 


2.25 


1.58 
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V. RECOMMENDATIONS 



GPS data from other days of the Phase II cruise should be used with the 
modified programs discussed here to further check our conclusions. 

Tests should be made in order to find the reason for offsets in a southeasterly 
direction and evaluate their correlation, if any, with the geometry of the satellites used. 

To test the conclusion that the accuracy depends on the large variations of the 
azimuth of the satellites, SVs other than SV6 and SVl 1 should be used. 

The relation between satellites' culminations during observations and the 
accuracy of the positions should be further investigated to provide a criterion for 
satellite selection. 

A better antenna site on the ship should be found to avoid shadowing and to 
minimize the effects of pitch and roll. Antennas of two similar systems should be 
installed at both high and low elevations to determine the effects of pitch and roll. 

The Kalman filter should be improved to avoid the oscillation when data are 
noisy. Other algorithms should be implemented and their results compared to the 
results of the KALMN2 program. 

To provide a real-time position computation the KALMX2 program or a similar 
program should be installed in a transportable computer. A system should be designed 
to have one computer converting and logging the data, another processing and 
displaying results, saving them or sending them to still another computer or data 
recording system, and there should be a master system to control the synchronization 
of time recorded with all data and to control the flow of information. 

Processing of the same data should be done using the precise ephemeris to see 
how the solution is improved. 

Processing in differential mode should be done to see how much the solution is 
improved. 

Comparison of data from other receivers should be done to check which type of 
receiver provides higher accuracy. 

Next an experiment should be designed to obtain point or differential positions in 
real time for moving platforms (survey ship) with little or no shore support. 
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To allow better evaluation of the GPS data, positions computed with more than 
two lines of position should be used as reference. 
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APPENDIX A 

PROGRAM FALCON. SOURCE LISTING 



//EZEQUIEL JOB ( 0812, 9999 EZECHJIEL CLASSIC 
//*MAIN 0RG=NPGVM1 . 0812P , LINES=( 99 ) ,CARDS=< 99 ) 

//■^^FORMAT PR,DDNAME=GO.FT06F001, 

//*F0RM$=SEP1 

// EXEC FORTVCLG 

//FORT.SYSIN DO * 

C 

C AUTHOR: AUGUSTO EZEQUIEL 

C DATE: APRIL 10, 1987 

C 

C DESCRIPTION: 

C THE PROGRAM TAKES THE RAW RANGE-RANGE DATA FROM AN INPUT FILE, 

C COMPUTS THE POSITIONS AND AND THE COURSE BETWEEN TWO ANTENNA 

C POSITIONS AND PRINTS THE UTM COORDINATES, THE COMPUTED COURSE 
C AND THE GPS TIME OF EACH OBSERVATION. 

C THE PROGRAM RUNS IN MVS 

C I/O SPECIFICATIONS: SEE END OF THIS JOB. 

C 

C THIS PROGRAM IS RESTRICT AS IS TO THE SEAFLOOR BENCHMARK EXPERIMENT 

C FOR OTHER USES DO NOT FORGET TO ENTER THE COORDINATES OF THE 

C STATIONS AND OTHER PARAMETERS. 

C 

C COORDINATES ARE IN WGS72 DATUM. 

C 

C THE PROGRAM WILL RUN WITH ANY AMOUNT OF DATA ONLY LIMITED TO DISK 
C SPACE 
C 

C ANY BLANK LINES WILL TERMINATE THE PROGRAM IN ERROR 

C 

C 

C 

C 

C 

REAL *8 X(2),Y(2),XPOS,YPOS,XOLD,YOLD 

INTEGER LADEG1,LADEG2,LAMIN1,LAMIN2,L0DEG1,L0DEG2,L0MIN1,L0MIN2 
INTEGER MONTH , YE AR , DAY , HOUR ,MIN , LADEG , LAMIN , LODEG , LOMIN 
REAL LASEC1,LASEC2,L0SEC1,L0SEC2,LASEG,L0SEG,SECS 
C 

C STATION 1 FERRIER 

C 

LADEG1=36 

LAMIN1=33 

LASEC1=53.748 

L0DEG1=121 

L0MIN1=S3 

L0SEC1=52.939 

C 

C STATION 2 DOME ECC 

C 

LA0EG2=36 

LAMIN2=18 

LASEC2=20.95 

L0DEG2=121 

L0MIN2=54 

LOSEC2=00.50 

C 

C COMPUTE X-Y FOR CONTROL STATIONS 
C 

CALL GPUTM(LADEG1,LAMIN1,LASEC1,LODEG1,LOMIN1,LOSEC1,Y(1) ,Xtl ) ) 
CALL GPUTM(LADEG2,LAMIN2,LASEC2,LODEG2,LOMIN2,LOSEC2,Y(2) ,X( 2 ) ) 
C 
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noo noo ooo o ooo ooo ooooo ooo oon ooo 



C ENG OF INITIAL DATA 

C 

C 

C COMPLTTS THE FIRST POSITION TO HAVE A STARTING POINT FOR 
C EVALUATION OF THE COURSE 

C 

READ(SO, * ,EN0=100) MONTH , GAY , YEAR , HOUR , MIN , SECS, 

KOGE 1 , RANGl , ISTRE 1 , KOGE 2 , RANG2 , ISTRE 2 
CALL RRXY( X,Y , RANGl ,RANG2,X0LG ,YOLG ) 

REAG IN GATA FROM INPUT FILE 

5G CONTINUE 

REAG(SG, * ,ENG=1G01 MONTH , GAY , YEAR , HOUR , MIN , SECS, 

* KOGE 1 , RANGl , ISTRE 1 , KOGE 2 , RANG2 , ISTRE 2 

COMPUTE X-Y ANG GP OF RR POSITION 

CALL RRXY( X,Y , RANGl ,RANG2 ,XPOS ,YPOS ) 

COMPUTES THE COURSE 

CALL HEA0( XOLG ,YOLG , XPOS,YPOS, COURSE , IE RR ) 

IF LANE JUMPS GOES NOT SAVE THE POSITION 

IF(IERR.LT.O) HRITE(6,1)IERR 
L FORMAT! 16) 

IF(IERR.LT.G) GOTO 50 

COMPUTES THE GEOGRAPHIC COOROINATES 

CALL UTMGP ( YPOS ,XPOS , LAOEG , LAMIN , LASEG , LOOEG , LOMIN , LOSEG ) 

OUTPLTT THE RESULTS 

WRITE! 6,60) MONTH , GAY , YEAR, HOUR, MIN, SECS, LAOEG, LAMIN, LASEG 
» , LOOEG , LOMIN , LOSEG , COURSE 

60 FORMAT! 2! IX, 12 ), IX, 14, 3! 1X,I3 ,1X,I2,1X,F7.3 ) ,1X,F5.0 ) 

WRITE ! 7 ) MONTH ,0AY , YEAR , HOUR ,MIN ,SECS , LAOEG , LAMIN , LASEG 
n , LOOEG, LOMIN, LOSEG 

SAVES THE COURSE IN A SEPARATE FILE 

WRITE ! 8 ) MONTH , GAY , YEAR , HOUR , MIN , SECS , COURSE 

SAVES POSITION FOR COURSE COMPUTATION 

XOLO=XPOS 
YOLO=YPOS 
GOTO 50 

ENG OF PROGRAM 

100 CONTINUE 
STOP 
ENG 
C 

SUBROUTINE RRXY!XR,YR,Rl,R2,XCO,YCO) 

c 

IMPLICIT REALMS !A-H,0-Z) 

OIMENSION XR!2),YR!2) 

A = OSQRT! !XR!2)-XR!l))9f*2 + ! YR! 2 )-YR! 1 ) )**2 ) 

CR = !XR! 2)-XR! 1) ) / A 
SR = !YR! 2)-YR! 1) ) / A 
XP = !R1*R1 - R2i^R2 + A^A ) / !2*A) 

ARG = Rl^Rl - XP^^XP 



66 



n o n o 



YP = OSQRT(ARG) 

XCO = XP*CR + YP*SR + XR(1) 

YCO = XP^SR - YP*CR + YR(1) 

RETURN 

END 

C 

^ \-/w\>w-\^ \^w\^\^wwv^wv/w-w«'\^w w v/wwwwwV^w wv#ww :^ wwwwwwwww 

SUBROUTINE GPUTM( LAOEG , LAMIN , LASEG , LOOEG , LOMIN , LOSEG , NORTH , EAST ) 
C 

DOUBLE PRECISION A ,R ,N ,AP ,BP ,CP ,0P ,EP ,S,R1 ,ESQ ,ESQP ,RM ,RP ,KO 
DOUBLE PRECISION R2,R3 ,R4,R5,P,P2,P3,P4,P5,P6 ,A6 ,B5,SINSEC 
DOUBLE PRECISION PHI ,DLAM, NORTH, EAST, B,PHIMIN, PI , LON 
REAL LASEG, LOSEG 

THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE ID CENTRAL MERIDIAN 123 DD DD W 



C 



C 

C 

C 



C 

C 

C 



c 



c 



c 



c 



c 



CM=-123.DDO 

PHI =DFLOAT( LADEG )+DFLOAT ( LAMIN )/6D . DD+DBLE ( LASEG )/36DD . DD 
LON=DFLOAT( LODEG )+DFLOAT( LOMIN )/6D . D0+DBLE( LOSEG )/360D . DD 
LON=-LON 

DLAM=( LON-CM 1^3600. ODD 

A=6378135.0D0 

R=298.26DD 

KO=D.9996DO 

B = A*(R-1.D0)/R 

N = (A-B)/(A+B) 

AP = A*( (1.00-N)+5.00/4.00*(N^2-N*^^3)+81.00/64.D0*(N**4-N**5) ) 

BP = 3 . DD/2 . D0*A*( ( N-N^^^2 ) + 7 . 00/8 . 00*( N**3-N**4 )+55 . 00/6^. 00^^N^^*5 ) 
CP = 15. DO/16. 00*A^^(N**2-N)f^^3+3. 00/4. 00*( N**4-N^t^^5 ) ) 

DP = 35. DO/48. 00«Ait(N**3-N**4+ll. 00/16. 00*N»*5) 

EP = 315. 00/512. 00i«A^^(N**4-N^t*5) 

PHIMIN = PHI*60. 00*2. 908882086660-4 

PI=0ARCOS(-1.00) 

PHI=PHI/180.00*PI 

S = AP*PHIMIN-BP*DSIN( 2.00*PHI)+CP*0SIN(4.00*PHI ) 

$ - 0P*0SIN( 6.00*PHI)+EP*DSIN(8.00*PHI ) 

R1 = KO*S 

SINSEC = ( 1.00/3600. 00 )/180.D0*PI 
SINSEC=OSIN(SINSEC) 

ESQ = ( A**2-B**2)/A**2 
ESQP = ESQ/( 1.00-ESQ) 

RM = A*( 1.00-ESQ )/(DSQRT( 1.00-ESQ*0SIN( PHI )**2) )**3 
RP = RM*(1.00+ESQP*DCOS(PHI)**2) 

R2 = RP*OSIN(PHI )*OCOS(PHI )*SINSEC**2/2.00*KO*1.008 
R3 = SINSEC**4*RP*0SIN(PHI)*0COS(PHI)**3/24.00*(5.00-0TAN(PHI)**2 
$ + 9.00*ESQP*0COS(PHI )**2+4.00*ESQP*ESQP*0COS(PHI )**4)*KO*1.016 
= RP*OCOS(PHI )*SINSEC*KO*1.04 
R5 = SINSEC**3*RP*DC0S( PHI )**3/6.00*( 1.00-0TAN( PHI )**2 
$ + ESQP*OCOS(PHI)**2 )*KO*1.012 

P = .000100*0LAM 
P2 = P**2 
P3 = P**3 

P4 = 

P5 = p**5 
P6 = P**6 

A6 = P6*SINSEC**6*RP*DSIN(PHI)*0COS(PHI)**5/720.D0 
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$ * ( 61. DO-58. D0*0TAN(PHI)**2+DTAN( PHI 

$ + 270 . DO*ESQP*DCOS( PHI )**2-330 . D0*ESQP*DSIN( PHI )**2 )*K0*1 .024 

B5 = P5*SINSEC**5*RP*DC0S( PHI )**5/120 . D0*( 5. DO-18. D0*DTAN( PHI )**2 
$ + DTAN(PHI)**4+14.DO*ESQP^^DCOS(PHI)**2 

$ - 58.D0*ESQP*0SIN(PHl)**2)*KO9^1.D20 

N0RTH=R1+R2*P2+R3*P4+A6 
EAST=(R4*P+RS*P3+B5)+500000.D0 
RETURN 
END 



SUBROUTINE HEAD( XOLD ,Y0LD ,XPOS , YPOS , COURSE , lERR ) 



10 

20 



REALMS XOLD , YOLD ,XP0S ,YPOS 
IERR=0 

DX=XP0S-X0LD 

DY=YP0S-Y0LD 

RADIUS=SQRT( DX**2+DY**2 ) 

IF( RADIUS. GT. 20. ) IERR=-1 
IF( RADIUS. GT. 0. ) GO TO 10 
COURSE =0.0 
RETURN 
CONTINUE 

COURSE = ARCOSCDY/RADIUS) 

COURSE =C0URSE*180 ./3 . 141592654 
IFCDX.LT.O.) COURSE =360. -COURSE 
RETURN 
END 



SUBROUTINE UTMGP ( NORTH , EAST , LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ) 



DOUBLE PRECISION A,R,N,B, BP ,CP ,DP ,EP ,S ,R1 ,ESQ,ESQP ,RM ,RP ,KO 
DOUBLE PRECISION R7 ,R8 ,R9 ,E5 ,Q ,Q2 ,Q3 ,Q4 ,Q5,Q6 ,D6 , AP ,SINSEC 
DOUBLE PRECISION RIO ,DLAM >NORTH , EAST ,RPHI ,DPHI ,PHIMIN 
DOUBLE PRECISION EPRIME ,DELTA,D6NUM,CM,PHI ,L0N,PI 
REAL LASEG, LOSEG 



THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN MGS 72 
IN ZONE 10 CENTRAL MERIDIAN 123 00 00 N 



CM=-123.0D0 



A=6378135.0D0 

R=298.26D0 



K0=0. 999600 

PI=DARCOS( -l.DO) 

B = A*(R-1.D0)/R 

N = (A-B)/(A+B) 

AP = A*( ( l.D0-N)+S.D0/4.D0*(N**2-N**3)+81.D0/64.D0*(N**4-N**5) ) 

BP = 3. DO/2. D0*A*( (N-N**2) + 7. DO/8. D0*(N**3-N**4)+55. DO/64. D0^^N**5) 
CP = 1S.D0/16.D0*A*(N**2-N**3+3.D0/4.D0*(N**4-N**5)) 

DP = 35. D0/48.D0*A*(N**3-N**4+11. DO/16. D0*N**5) 

EP = 315.D0/512.D0*A*(N**4-N**5) 

FIRST APROXIMATION OF PHI 

PHI =NORTH/30 . 8D0/3600 . DO 

COMPUT TRUE MERIDIONAL DIST AND APROXIMATE PHI 



DO 100 1=1,8 

PHIMIN=PHI*60 . D0*2 . 90888208666D-4 
RPHI=PHI/180.D0*PI 
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S = APi^PHIMIN-BP*DSIN( 2.D0*RPHI )+CP*DSIN(4.D0*RPHI ) 
$ -DP*DSIN(6.D0*RPHI )+EP*DSIN( 8 . DO^^RPHI ) 

R1 = KO*S 
DELTA=N0RTH-R1 
PHI=( DELTA/20.8/3600. )+PHI 
100 CONTINUE 

RPHI=PHI/180.D0*PI 



SINSEC = ( l.D0/3600D0)/180.D0*PI 
SINSEC = DSIN( SINSEC) 

C 

ESQ = (A^h^2-B**2)/A**2 
ESQP = ESQ/( 1. DO-ESQ) 

RM = A*(l. DO-ESQ )/(DSQRT(l.D0-ESQ*DSIN(RPHI)**2))**3 
C 

RP = RM*(1.DO+ESQP*DCOS(RPHI)^«^2) 

R7=( DTANl RPHI )/( 2 . DO*RP*9^2^^SINSEC ) )*( 1 . DO+ESQP*DCOS( RPHI )^Z ) 

$ *l.D12/KO**2 

R8 = ( DTAN( RPHI )/( 2<+. D0*RP^e*4*SINSEC ) )*( 5 . DO+3 . DO*DTAN( RPHI )**2 
$ +6 . DO*ESQP*DCOS( RPHI )**2-6 . DO*ESQP*DSIN( RPHI )^«^2-3 . D0*ESQP**2 

$ *DCOS( RPHI )i^*4-9.D0*ESQP**2*DC0S( RPHI )^«t2*DSIN( RPHI )**2 ) 

$ *1.D24/K0**4 

R9=l . D0/DC0S( RPHI )/( RP*SINSEC )*1 . D6/K0 

R10 = l . DO/DCOS( RPHI )/( 6 . D0*RP**3*SINSEC )*( 1 . DO+2 . DO*DTAN( RPHI 
$ ♦ESQP?^DCOS(RPHI )**2 )^^1 . D18/K0**3 

C 

EPRIME=EAST-500000 . DO 

Q = .000001D0*EPRIME 

Q2 = Q**2 

Q3 = Q**3 

Q4 = Q**4 

Q5 = Q**5 

Q6 = Q**6 

D6NUM=Q6*DTAN( RPHI ) 

D6 = ( D6NUM )/( 720.D0*RP**6*SINSEC )*( 61 . DO + 90 . D0*DTAN( RPHI )**2 
$ +45 . DO*DTAN( RPHI )i^*4 + 107 . D0*ESQP*DC0S( RPHI )**2 

$ -162 . DO*ESQP*DSIN( RPHI )**2-45 . DO*ESQP*DTAN( RPHI )**2 

$ *DSIN( RPHI )^Z )*1 . D36/K0**6 
E5=( Q5*l . DO/DCOS( RPHI ) )/( 120 . D0*RP**5*SINSEC W 5 . DO + 28 . DO 
$ *DTAN( RPHI )**2 + 24 . DO*DTAN( RPHI )**4+6 . DO*ESQP*DCOS( RPHI )¥^Z 

$ +8 . DO*ESQP*DSIN( RPHI )**2 )*1 . D30/K0**5 

C 

DPHI = ( -R7i^Q2+R8*Q4-D6 )/3600 . DO 
DLAM=( R9*Q-R10*Q3+E5 )/3600. DO 
PHI=PHI+DPHI 
LON=CM+DLAM 

CALL OMS(PHI,LADEG,LAMIN,LASEG) 

CALL DMS( LON,LODEG,LOMIN,LOSEG) 

RETURN 

END 

C 



SUBROUTINE DMS( DEC,LDEG,MIN,SEC ) 



C 

C 



DOUBLE PRECISION DEC ,XNUM,XMIN 



XNUM=DABS( DEC ) 

LDEG=DINT(XNUM) 

XMIN=(XNUM-DFLOAT( LDEG) )*60.D0 
MIN=DINT(XMIN) 

XNUM=( XMIN-DFLOAT( MIN ) )*60 . DO 
SEC=SNGL(XNUM) 

IF (LDEG.GE.360) LDEG=LDEG-360 

RETURN 

END 



/J^ 

//GO.FT06F001 DD SYSOUT=* 

//GO.FT07F001 DD DSN=MSS . S0812 . FALCON. ANTENNA . POS ,DISP=OLD 
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//GO.FT08F001 DD 
//GO.FT50F001 DD 



// 


DO 


// 


DD 


// 


DO 


// 


DD 


// 


DO 


// 


DD 


// 





DSN=MSS . S0812 . FALCON . COURSE . ROUGH ,DISP=OLO 
OSN=MSS. S081Z . FALCON . OATA( FALCONl ) ,DISP=SHR 
DSN=MSS.S0812.FALC0N.0ATA( FALCON2 ) ,DISP=SHR 
DSN=MSS. S0812 . FALCON . DATA( FALCONS ) ,DISP=SHR 
DSN=MSS.S0812 . FALCON. DATA( FALC0N4 ) ,DISP=SHR 
0SN=MSS.S0812. FALCON. DATAt FALCONS ) ,DISP=SHR 
DSN=MSS.S0812. FALCON. DATA( FALC0N6 ) ,DISP=SHR 
DSNsMSS . S0812 . FALCON . DATA( FALCON7 ) ,0ISP=SHR 
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APPENDIX B 

PROGRAM COURSE SMOOTH. SOURCE LISTING 



//JOBSMOTH JOB ( 0812,9999 ), ’EZEQUIEL ' ,CLA$S=C 
//*MAIN ORG=NPGVMl . 0812P ,LINES=( 99 ) ,CARDS=( 99 ) 

//^^FORMAT PR,DDNAME=GO.FT06F001, 

//*FORMS=SEPl 

// EXEC FORTVCLG 

//FORT.SYSIN DD * 

C 

C AUTHOR: AUGUSTO EZEQUIEL 
C DATE; JANUARY 05, 1986 
C 

C DESCRIPTION: 

C THE PROGRAM MAKES THE RUNNING AVERAGE OF THE COURSE COMPUTED 

C BY THE PROGRAM < COMP > IN ORDER TO TAKE OUT THE EFECTS OF THE 

C PITCH AND ROLL EFFECTS. THE FIRST AND LAST TWO POSITIONS ARE 
C REJECTED AS NO RUNNING AVERAGE IS POSSIBLE FOR THEM. 

C THE PROGRAM RUNS IN MVS 

C I/O SPECIFICATIONS: SEE END OF THIS JOB. 

C 

C THIS PROGRAM IS RESTRICT AS IS TO THE SEAFLOOR BENCHMARK EXPERIMENT 



C 

C THE PROGRAM WILL RUN WITH ANY AMOUNT OF DATA ONLY LIMITED TO DISK 
C SPACE 
C 

C ANY BLANK LINES WILL TERMINATE THE PROGRAM IN ERROR 

C 

C 

C 

C 

C 

REAL C0URSE(5),SECS(5) 

INTEGER MONTH( 5 ) ,YEAR( 5 ) ,DAY( 5 ) ,HOUR( 5 ) ,MIN( 5 ) ,FLAG 
LOGICAL SAVE 
C 

C INITIALIZATION OF FLAGS 

C 

SAVE= .FALSE. 

FLAG=1 

C 

C READS THE FIRST FIVE COURSES 

C 

DO 100 J=l,5 

READ( 51 ,END=1100 ) MONTH( J ) ,DAY( J ) ,YEAR( J ) ,HOUR( J ) , 

$ MIN( J),SECS( J),COURSE( J) 

150 FORMAT( 2( IX ,12 ) , IX ,14 ,1X ,12 , 1X,I2 , 1X,F4. 1 , IX, F07 . 2 ) 

C WRITE(6,150) MONTH(J),DAY( J),YEAR( J),HOUR( J), 

C $ MIN( J),SECS(J ), COURSE! J) 

100 CONTINUE 
C 

C AVERAGE THE COURSE 

C 

105 AVG=( COURSE! 1 l+COURSE! 2 )+COURSE! 3 l+COURSE ! 4 )+COURSE! 5 ) )/5. 0 
C 

C SAVES THE AVERAGE IN THE MIDDLE POSITION 

C 

COURSE! 3 )=AVG 



C 

C 

C 

C 

C 



DOES NO SAVE THE FIRST TWO POSITONS OF THE FILE AS NO RUhWING 
AVERAGE EXISTS FOR THEM 

IF! SAVE) GOTO 50 
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FLAG=FLAG+1 

IF( FLAG.EQ.3) SAVE=.TRUE. 
GOTO 70 



C 

C WRITES THE DATA WITH THE AVERAGED COURSE INTO THE FILE 
C 

50 WRITE( 7 ) MONTHf 1 ) ,DAY( 1 ) >YEAR( 1 ) >HOUR( 1 ) ,MIN( 1 ) ,SECS( 1 ) , COURSE ( 1 ) 
C 

WRITER 6 ,150 ) MONTHS 1 ) ,DAY( 1 ) ,YEAR( 1 ) ,HOUR( 1 ) , 

$ MIN( 1 ) ,SECS( 1 ) ,COURSE( 1 ) 

C 

C MOVES FORWARD THE DATA 
C 

70 CONTINUE 

DO 75 J=2,5 
MONTH! J-1)=M0NTH(J) 

DAY! J-1)=DAY(J) 

YEAR! J-1)=YEAR( J) 

HOUR! J-1)=H0UR! J ) 

MIN! J-1)=MIN( J) 

SECS! J-1)=SECS! J) 

COURSE! J-1 )=COURSE! J) 

75 CONTINUE 
C 

C READS ONE MORE COURSE 
C 

READ! 51 ,END=1000 ) MONTH! 5 ) ,DAY! 5 ) ,YEAR! 5 ) ,HOUR! 5 ) , 

$ MIN! 5), SECS! 5), COURSE! 5) 

C 

C GO BACK TO COMPUT A NEW AVERAGE 
GO TO 105 
C 

C END OF FILE FOUND DURING THE INPUT OF A FITH COURSE 
C SO SAVES THE REMAINING TWO AVERAGES 

1000 CONTINUE 

DO 1001 J=l,2 

WRITE! 7,150 ) MONTH! J ) ,DAY! J ) ,YEAR! J ) ,HOUR! J ) , 

$ MIN! J), SECS! J), COURSE! J) 

WRITE ! 6 , 150 ) MONTH! J ) , DAY! J ) ,YEAR! J ) ,HOUR! J ) , 

$ MIN! J), SECS! J), COURSE! J) 

1001 CONTINUE 

STOP 



C 

C NO OUTPUT IS MADE IF AN END OF FILE WAS FOUND DURING THE READING 

C OF THE FIRST FIVE COURSES 

1100 CONTINUE 

WRITE!6,1<^0) 

140 FORMAT !* ERROR. END OF FILE DURING THE FIRST SET OF COURSES') 
STOP 
END 

/* 



//G0.FT06F001 

//G0.FT07F001 

//GO.FT51F001 

// 



DD SYSOUT=« 

DD DSN=MSS . S0812 . FALCON. COURSE . SMOOTH ,DISP=SHR 
DD DSN=MSS . S0812 . FALCON . COURSE . ROUGH , DISP=SHR 
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APPENDIX C 

PROGRAM TRANSDUC FALCON. SOURCE LISTING 



//JOBTRANS JOB ( 0812, 9999 ),' EZEQUIEL ‘ ,CLASS=C 

//^MAIN 0RG=NPGVM1.0812P,LINES=( 99),CARDS=( 99) 

//^FORMAT PR,0DNAME=G0.FT06F001, 

//*F0RMS=SEP1 

// EXEC FORTVCLG 

//FORT.SYSIN DD * 

C 

C 

C AUTHOR: AUGUSTO E2EQUIEL 

C DATE; MARCH 04, 1987 

C 

C DESCRIPTION: 

C 

C THIS PROGRAMS TAKES EACH POSITION OF THE ANTENNA 

C THE COMPUTED COURSE ( FILTERED OR NOT ) ,THE PITCH 

C AND ROLL DATA, AND COMPUTES THE POSITIONS OF THE 

C TRANSDUCER, APLYING THE SEVEN PARAMETER TRANSFORMATION 

C TO THE OFFSETS OF THE TRANSDUCER IN RELATION TO THE ANTENNA 

C CONSIDERED THE CENTER OF COORDINATE SYSTEM AND THEN ADDING 

C THE CORRECTIONS TO THE COORDINATES. 

C 

C THE PROGRAM RUNS IN MVS 

C 

C I/O SPECIFICATIONS: SEE END OF THIS JOB. 

C 

C THIS PROGRAM IS RESTRICT AS IS TO THE SEAFLOOR BENCHMARK 
C EXPERIMENT 

C 

C THE PROGRAM WILL RUN WITH ANY AMOUNT OF DATA ONLY LIMITED TO DISK 
C SPACE 

C 

C ANY BLANK LINES WILL TERMINATE THE PROGRAM IN ERROR 
C 
C 
C 

DOUBLE PRECISION XP0S,YP0S,TIME1 ,TIME2,TIME3,SECD,SECH,RATE 
DOUBLE PRECISION TIME4, TIMES 

REAL COURSE ,PITCH( 2 ) ,ROLL( 2 ) ,SEC1( 2 ) ,SEC2,PIT ,ROL ,SEC3( 2 ) 

REAL DX,DY,DZ,0FFX,0FFY,0FF2,HEAD( 2),LASEG,L0SEG 
INTEGER MONTHl ( 2 ) ,DAY1 ( 2 ) ,YEAR1( 2 ) ,H0UR1 ( 2 ) ,MIN1 ( 2 ) 

INTEGER M0NTH2 ,DAY2 ,YEAR2 ,H0UR2 ,MIN2 
INTEGER M0NTH3( 2),DAY3( 2),YEAR3( 2),H0UR3( 2),MIN3( 2) 

INTEGER LADEG,LAMIN,L0DEG,L0MIN 
C 

C INITIALIZATION OF CONSTANTS 
C 

SECH=3600.D0 

SECD=24.D0*SECH 

C 

0FFX=-5.34S 

0FFY=-9.4S3 

0FFZ=-16.1S2 

C 

C READS TWO SETS OF PITCH AND ROLL DATA 
C 

READ(S1,*,END = 100) MONTHK 1 ) ,DAY1 ( 1 ) ,YEAR1 ( 1 ) ,H0UR1( 1) ,MIN1( 1 ) , 
$ SECK 1 ) ,PITCH( 1 ) ,R0LL( 1 ) 

TIME1=DFL0AT( DAYK 1 ) )*SECD+DFLOAT( HOURl ( 1 ) )*SECH 
$ + DFLOAT(MIN1(1))*60.DO+DBLE(SEC1(1)) 

C 

READ(S1,*,END = 100) MONTHl ( 2 ) ,DAY1( 2 ) , YEARl ( 2 ) ,H0UR1( 2 ) ,MINl( 2 ), 



73 



ooo ooo ooo ooo oooooo o ooo o non 



$ SECK 2 ), PITCH! 2 ) ,ROLU 2 ) 

TIME2=DFL0AT( DAYK 2 ) )J^SECD+OFLOAT( HOURK 2 ) )J^SECH 
$ + DFLOAT(MIN1(2))^^60.DO+OBLE(SEC1(2)) 

READS TWO SETS OF COURSE 

READ! 52 ,END=100 ) MONTHS! 1 ) ,0AY3! 1 ) , YEARS! 1 ) , HOURS! 1 ) ,MINS! 1 ) > 

$ SECS! 1 ), HEAD! 1) 

TIME4=DFL0AT! DAYS! 1 ) IJ^SECD+DFLOAT! HOURS! 1 ) )*SECH 
$ + DFLOAT!MINS!l))i^60.D0+DBLE!SECS! 1)) 

READ! 52 ,END = 100 ) MONTHS! 2 ) ,DAYS! 2 ) , YEARS! 2 ) >HOURS! 2 ) ,MINS! 2 ) , 

$ SECS! 2), HEAD! 2) 

TIMES=DFLOAT! DAYS! 2 ) ))^SECD+DFLOAT! HOURS! 2 ) )^^SECH 
$ + DFLOAT!MINS!2) )i^60.D0+DBLE!SECS!2)) 

READS ONE POSITION OF THE SHIP 

20 READ!50,END=100) MONTH2 >DAY2 >YEAR2 ,HOUR2 ,MIN2 ,SEC2 , 

$ LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG 

TIMES=DFLOAT! DAY2 )#SECD+DFLOAT! HOUR2 )^^SECH 
$ + DFLOAT!MIN2)*60.D0+DBLE!SEC2) 

POSITION IN TIME BETWEEN THE TWO SETS OF PITCH AND ROLL DATA 
AND THE TWO SETS OF COURSE DATA 
25 CONTINUE 

POSITION TO EARLY FOR PITCH AND ROLL DATA 
IF!TIMES.GE. TIMED GO TO 40 

WRITE ! 6 ,S0 ) MONTH2 >DAY2 ,YEAR2 >H0UR2 ,MIN2 ,SEC2 
SO FORMAT!' PSN 3 ’ ,12 ,1X,I2 ,1X,I4 ,1X,I2 ,1X,I2 ,1X,F4. 1 , 

$ ' REJECTED. NO PITCH AND ROLL DATA') 

GOTO 20 

POSITION TO EARLY FOR PITCH AND ROLL DATA 
40 CONTINUE 

IF!TIMES.GE.TIME4) GO TO 60 

WRITE! 6,50 ) MONTH2,DAY2 ,YEAR2,HOUR2,MIN2,SEC2 
so FORMAT!' PSN 3 ' ,12 ,1X ,12 , 1X,I4 ,1X ,12 , IX ,I 2 , IX, F4. 1 , 

$ • REJECTED. NO COURSE DATA') 

GOTO 20 

POSITION LATER THAN THE SECOND SET OF PITCH AND ROLL DATA 
60 CONTINUE 

IF!TIMES.LE.TIME2) GOTO 70 

MOVES THE PITCH AND ROLL DATA TO THE FIRST SET 

M0NTH1!1)=M0NTH1!2) 

DAY1!1)=DAY1!2) 

YEAR1!1)=YEAR1!2) 

H0UR1!1)=H0UR1! 2) 

MIN1!1)=MIN1!2) 

SECl! 2)=SEC1! 2) 

PITCH!1)=PITCH! 2) 

ROLL!l)=ROLL!2) 

TIME1=TIME2 

READS A SECOND SET OF PITCH AND ROLL DATA 

READ! 51 ,END = 100 ) MONTHl! 2 ) ,DAY1! 2 ) ,YEAR1! 2 ) ,H0UR1! 2 ) ,MIN1! 2 ) , 
$ SECl! 2), PITCH! 2), ROLL! 2) 

TIME2=DFLOAT! DAYl! 2 ) )*SECD + DFLOAT! HOURl! 2 ) )*SECH 
$ + DFLOAT!MIN1!2))^^60.DO+DBLE!SEC1!2)) 

GOTO 25 
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POSITION LATER THAN THE SECOND SET OF COURSE DATA 
70 CONTINUE 

IF( TIMES. LE. TIMES) GOTO 80 

MOVES THE COURSE DATA TO THE FIRST SET 

MONTHS! l)=MONTH3( 2) 

DAY3(1)=DAY3(2) 

YEAR3(1)=YEAR3(2) 

HOURS! 1 )=HOUR3! 2 ) 

MINS! 1)=MIN3!2) 

SEC3!2)=SEC3! 2) 

HEAD!1)=HEAD!2) 

TIME4=TIME5 

READS A SECOND SET OF PITCH AND ROLL DATA 

READ! 52 ,END = 100 ) MONTHS! 2 ) ,DAY3! 2 ) , YEARS! 2 ) , HOURS! 2 ) ,MIN3! 2 ) , 

$ SECS! 2), HEAD! 2) 

TIMES=DFLOAT! DAYS! 2 ) )*SECD+0FLOAT! HOURS! 2 ) )*SECH 
$ + DFLOAT!MIN3!2))*60.D0+DBLE!SEC3!2)) 

GOTO 2S 

POSITION HITHIN THE TWO SETS IN TIME 

80 RATE=! TIME3-TIME1 )/! TIME2-TIME1 ) 

PIT=PITCH! 1 ) + ! PITCH! 2 )-PITCH! 1 ) )*RATE 
ROL=ROLL! 1 ) + ! ROLL! 2 )-ROLL! 1 ) )^RATE 
RATE=! TIME3-TIME4 )/! TIMES-TIME^ ) 

COURSE = HEAD! 1 ) + ! HEAD! 2 )-HEAD! 1 ) )^fRATE 

COMPUTS THE CORRECTIONS AND THE POSITIONS OF THE TRANSDUCER 

CALL PTH!OFFX,OFFY ,OFFZ, COURSE ,PIT ,ROL ,DX,DY ,DZ ) 

COMPUTES THE UTM COORDINATES 

CALL GPUTM! LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG , YPOS ,XPOS ) 

XPOS=XPOS+DX 

YPOS=YPOS+DY 

CONVERTS BACK TO GP 

CALL UTMGP ! YPOS ,XPOS , LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ) 
OUTPUTS THE RESULT 

WRITE! 7 ) MONTH2 ,DAY2 ,YEAR2 ,HOUR2 ,MIN2 ,SEC2 , 

LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG 

WRITE ! 6 , 90 ) MONTH2 , DAY2 ,YEAR2 ,HOUR2 ,MIN2 ,SEC2 , 
n LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG 

90 FORMAT! 2! 1X,I2 ),1X,I<^,1X,I2,1X,I2,1X,F^^. 1,2! 1X,I3,1X,I2,1X,F7.4) ) 

READS ONE MORE POSITION 
GOTO 20 

END OF PROGRAM 

100 CONTINUE 
STOP 
END 

SUBROUTINE PTH! XOFF ,YOFF ,20FF , COURSE , PITCH , ROLL ,DX,DY ,DZ ) 

THIS SUBROUTINE COMPUTES THE CORRECTIONS TO THE 
COORDINATES DUE TO THE OFSET OF ANTENA IN RELATION 
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TO THE TRANSDUCER 

PI=ARCOS(-1.0) 

ROLL=ROLL^^PI/180. 

PITCH=PITCH*PI/180. 

COURSE=COURSE*PI/180 . 

CHANGE SIGNAL OF ORIGINAL DATA TO MATCH THE RIGHT HAND 
SYSTEM CONVECTION 

ROLL=-ROLL 

COSPHI= COS(ROLL) 

SINPHI= SIN(ROLL) 

COSK= COS (COURSE) 

SINK= SIN( COURSE) 

COSN= COS( PITCH) 

SINNs SIN( PITCH) 



ox 


= XOFF 




COSPHI 


* COSK + 












$ 


YOFF 




( COSW 


* SINK + 


SINW 


* SINPHI 


* 


COSK 


) + 


$ 


ZOFF 


* 


( SINW 


* SINK - 


COSW 


* SINPHI 


* 


COSK 


) 


DY 


= XOFF 


* 


COSPHI 


^ SINK * 


(-1. 


) + 








$ 


YOFF 


* 


( COSW 


* COSK - 


SINW 


* SINPHI 


* 


SINK 


) + 


$ 


ZOFF 


* 


( SINW 


* COSK + 


COSW 


* SINPHI 




SINK 


) 


DZ 


= XOFF 


* 


SINPHI 














$ 


YOFF 


* 


( SINW * COSPHI 


: ) + 










$ 


ZOFF 




( COSW 


* COSPHI 


) 











RETURN 

END 



SUBROUTINE UTMGP ( NORTH , EAST , LAOEG > LAMIN , LASEG > LODEG , LOMIN , LOSEG ) 

DOUBLE PRECISION A,R,N,B, BP ,CP ,DP ,EP ,S ,R1 ,ESQ,ESQP ,RM ,RP ,K0 
DOUBLE PRECISION R7 ,R8 ,R9 ,E5 ,Q ,Q2 ,Q3 ,Q4 >(?5 >(?6 ,D6 >AP ,SINSEC 
DOUBLE PRECISION RIO ,DLAM, NORTH , EAST ,RPHI ,DPHI ,PHIMIN 
DOUBLE PRECISION EPRIME ,DELTA,06NUM,CM,PHI ,LON,PI 
REAL LASEG, LOSEG 

THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 123 00 00 W 

CM=-123.000 

A=6378135.000 

R=298.26D0 

KO=0. 999600 

PI=0ARCOS(-1.00) 

B = A*(R-1.00)/R 

N = (A-B)/(A+B) 

AP = A*( (l.OO-N)+5.DO/‘^.O0*(N**2-N**3)+81.DO/64.OO*(N**4-N**5)) 

BP = 3.00/2.00*A*( (N-N**2) + 7.00/8.D0*(N**3-N**4)+S5.00/6‘^.D0*N**5) 
CP = IS. DO/16. 00*A*(N**2-N**3+3. DO/4. 00*( N**4-N**5 ) ) 

OP = 35. 00/48. 00*A*(N**3-N**4+11. 00/16. D0*N**5) 

EP = 315. 00/512. D0*A*(N**4-N**5) 

FIRST APROXIMATXON OF PHI 

PHI =NORTH/30 . 800/3600 . 00 

COMPUT TRUE MERIDIONAL OIST AND APROXIMATE PHI 



00 100 1 = 1,8 
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PHIMIN=PHI3^60 . D0*2 . 90888208666D-<+ 

RPHI=PHI/180 . DO*PI 

S = AP*PHIMIN-BP*DSIN( 2.D0*RPHI)+CP*DSIN(^.D0*RPHI) 
$ -DP*0SIN( 6.D0*RPHI )+EP#DSIN( 8 . 00*RPHI ) 

R1 = KO*S 
DELTA=N0RTH-R1 
PHI=(0ELT A/30. 8/3600. )+PHI 
100 CONTINUE 

RPHI=PHI/180.D0*PI 



SINSEC = ( 1.D0/3600D0 )/180.D0*PI 
SINSEC = 0SIN( SINSEC) 

ESQ = ( A**2-B**2 )/A^«(2 
ESQP = ESQ/( 1. DO-ESQ) 

RM = A*( 1, DO-ESQ )/( DSQRK l.D0-ESQ*DSINt RPHI )**2) )**3 
RP = RM#(1.DO+ESQP*DCOS(RPHI)**2) 

R7=( DT AN( RPHI )/( 2 . D0*RP**2*SINSEC ))*(!. D0+ESQP*DC0S( RPHI ) 

$ *l.D12/KO**2 

R8=( DTAN( RPHI )/( 2^ . D0*RP**4*SINSEC ) )*( 5. DO+3 . D0*DTAN( RPHI )**2 
$ +6 . DO*ESQP*DCOS( RPHI )**2-6 . DO*ESQP*DSIN( RPHI )**2-3 . D0*ESQP**2 

$ J^DCOS( RPHI )**^-9.DO*ESQP**2*DCOS( RPHI )**2*DSIN( RPHI )**2 ) 

$ *1.D24/K0**4 

R9=l . DO/DCOS( RPHI )/( RP*SINSEC )*1 . D6/KO 

R10=l . DO/DCOS( RPHI )/( 6 . DO*RP**3*SIf4SEC )*( 1 . DO+2 . D0*DTAN( RPHI )#*2 
$ +ESQP*DCOS( RPHI )^Z )*1 . D18/KO**3 

EPRIME=EAST-500000 . DO 

Q = .000001D0*EPRIME 

Q2 = Q**2 

Q3 = Q**3 

Q<+ = Q**4 

QS = Q5«^5 

Q6 = Q**6 

D6NUM=Q6#DTAN( RPHI ) 

D6 = ( D6NUM)/( 720,D0*RP**6*SINSEC )*( 61 . DO + 90 . D0*DTAN( RPHI )^Z 
$ +<+5.D0*DTAN(RPHI )**<++107. DO*ESQP*DCOSt RPHI )**2 
$ -162.D0*ESQP*DSIN( RPHI )**2-^5 . DO*ESQP*DTAN( RPHI )**2 
$ J^DSIN( RPHI )**2 )*1 . D36/KO**6 
ES = ( QS^l . DO/DCOS( RPHI ) )/( 120 . D0*RP**5*SINSEC )*( 5 . DO + 28 . DO 
$ *DTAN( RPHI )*#2+2<+. DO*DTAN( RPHI )^Ht4+6 . DO*ESQP*DCOS( RPHI )^Z 
$ +8 . DO*ESQP*DSIN( RPHI )**2 )*1 . OZO/KO^S 

DPHI=( -R7*Q2+R8*Q4-D6 )/3600.D0 
DLAM=( R9*Q-R10*Q3+E5)/3600.D0 
PHI=PHI+DPHI 
LON=CM+DLAM 

CALL DMS(PHI,LADEG,LAMIN,LASEG) 

CALL DMS( LON,LODEGaOMIN,LOSEG) 

RETURN 

END 

SUBROUTINE GPUTM( LADEG > LAMIN , LASEG , LODEG > LOMIN , LOSEG .NORTH . E AST ) 

DOUBLE PRECISION A ,R ,N >AP ,BP ,CP ,DP ,EP ,S,R1 ,ESQ,ESQP ,RM,RP ,KO 
DOUBLE PRECISION R2 ,R3 ,R<+,R5 ,P ,P2 ,P3 ,P<+, PS ,P6 .A6 >B5 , SINSEC 
DOUBLE PRECISION PHI .DLAM^NORTH .EAST ,B ,PHIMIN ,PI , LON 
REAL LASEG .LOSEG 

THIS SUBROUTINE CONPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 123 00 00 H 

CM=-123.0D0 

PHI=DFLOAT( LADEG ) + DFLOAT( LAMIN )/60 .D0+DBLE( LASEG )/3600 . DO 
LON=DFLOAT( LODEG )+DFLOAT( LOMIN )/60 . D0+D8LE( LOSEG )/3600 . DO 
LON =- LON 

DLAM=( LON-CM)*3600.0D0 
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A=6378155.000 

R=298.26D0 

C 

KO=0.9996D0 



C 



8 = A*(R-1.D0)/R 



C 



C 

C 

C 



C 

c 



c 



c 



c 



c 

c 

c 






N = (A-B)/(A+B) 

AP = A*( ( 1.00-N)+5.00/^♦.00*(N**2-N**3)+81.00/6<^.00*(N**4-N**5n 
BP = 3.00/2.D0*A*( (N-N**2 )+7. 00/8. 00*( N**3-N**^ )+55. 00/6^. D0*N**5 ) 
CP = 15.00/16.00*A*(N**2-N**3+3.00/^.00*(N**^-N**51) 

OP = 35. 00/48.00*A*(N**3-N**4+ll. 00/16. D0*N**5) 

EP = 315. 00/512. 00*A*(N**^-N**5) 

PHIMIN = PHI*60. 00*2. 908882086660-4 

PI=OARCOS( -1.00) 

PHI=PHI/180.00*PI 

S = AP*PHIMIN-BP*OSIN( 2.00*PHI )+CP*0SIN(4.00*PHI 1 
$ - DP*DSIN(6.00*PHI)+EP*0SIN(8.D0*PHI) 

R1 = KO*S 

SINSEC = C 1.00/3600.00 1/180. 00*PI 
SINSEC=0SIN(SINSEC1 

ESQ = ( A**2-B**2 )/A**2 
ESQP = ESQ/( 1.00-ESQ 1 

RM = A*( 1.00-ESQ 1/(0SQRTC1.00-ESQ*0SIN( PHI 1**2 )1**3 
RP = RM*( 1.00+ESQP*DC0S( PHI 1**21 

R2 = RP*0SIN(PHI 1*DCOS(PHI1*SINSEC**2/2.00*KO*1.008 
R3 = SINSEC**4*RP*DSIN(PHI1*OCOS(PHI1**3/24.00*(5.00-OTAN(PHI1**2 
$ + 9.DO*ESQP*DCOS(PHI1**2+4.DO*ESQP*ESQP*OCOS(PHI1**41*KO*1.016 
= RP*DCCS(PHI 1*SINSEC*K0*1.04 
R5 = SINSEC**3*RP*0C0S( PHI 1**3/6.00*(1.00-0TAN( PHI 1**2 
$ + ESQP*OCOS( PHI 1**2 1*K0*1. 012 

P = .000100*0LAM 
P2 = P**2 
P3 = P**3 

P4 = 

P5 = P**5 
P6 = P**6 

A6 = P6*SINSEC**6*RP*0SIN( PHI 1*0C0S( PHI 1**5/720.00 
$ * (61.00-58.00*0TAN( PHI 1**2+0TAN( PHI 1**4 

$ + 270 . 00*ESQP*0C0S ( PHI 1**2-330 . 00*ESQP*0SIN( PHI 1**2 1*K0*1 . 024 

B5 = P5*SINSEC**5*RP*0C0S( PHI 1**5/120.00*(5,00-18.00*OTAN( PHI 1**2 
$ + DTAN( PHI 1**4+14.00*ESQP*0C0S( PHI 1**2 

$ 58.00*ESQP*0SIN( PHI 1**2 1*K0*1. 020 

N0RTH=R1+R2*P2+R3*P4+A6 
EAST=( R4*P+R5*P3+B5 1+500000 . DO 
RETURN 
ENO 

SUBRCXJTINE OMS( OEC , LOEG ,MIN ,SEC 1 

OOUBLE PRECISION DEC ,XNUM,XMIN 

XNUM=DABS(0EC1 

L0EG=0INT(XNUM1 

XMIN=(XNUM-OFLOAT( LOEG 1 1*60. DO 
MIN=0INT(XMIN1 

XNUM=( XMIN-DFLOAT ( MIN 1 1*60 . 00 
SEC=SNGL(XNUM1 

IF (LDEG.GE.3601 L0EG=L0EG-360 

RETURN 

ENO 
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//GO.FT06F001 DO 
//GO.FT07F001 DO 
//GO.FT50F001 DD 
//GO.FT51F001 DD 
//GO.FT52F001 DD 
// 



SYSOUT=J^ 

DSN=MSS . S0812 . FALCON . TRANSDUC . POS ,DISP=SHR 
DSN=MSS . S0812 . FALCON . ANTENNA . POS >DISP=SHR 
DSN=MSS . S0812 . PITROLL ,DISP=SHR 
DSN=MSS . S0812 . SDAS . COURSE , DISP=SHR 
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APPENDIX D 

PROGRAM PLOT. SOURCE LISTING 

//EZEQUIEL 
// 

//FORT.SYSIN 
C 

C PROGRAM 

C 

C RUNS IN 

C 

C AUTHOR: 

C 

C DATE : 26 JANUARY 1987 

C 

C THIS PROGRAM MAKES A PLOT OF THE POSITIONS OF THE SHIP 

C 

C 

DOUBLE PRECISION XPOS,YPOS,XLEFT ,YLEFT,TIME1,TIME2,SECD,SECH 
REAL SECS, XPLT,YPLT,BLXH,BLYH, SCALE, VALUE, LASEG,LOSEG 
INTEGER MONTH , YEAR , DAY , HOUR , MIN , IPEN , LADEG , LAMIN , LODEG , LOMIN 

PLOTTER INITIALIZATION 

CALL PLOTS (0,0,0) 

CALL PLOT (2. ,2., -3) 

IPEN=3 

DIMENSIONS OF SHEET, LEFT CORNER AND SCALES 

BLXH=36. 

BLYH=^9. 

XLEFT=564800.D0 
YLEFT=4037000.D0 
SCALE=1./10000. 

INITIALIZATION OF CONSTANTS 

TIME1=0.D0 
SECH=3600.D0 
SECD = 2<+.D0*SECH 

GRID 



CALL 


PLOT 


(0.1, 0.1, 3) 




CALL 


PLOT 


(0.1,BLYH-0.1, 


2) 


CALL 


PLOT 


(BLXH-0.1,BLYH 


-0.1,2) 


CALL 


PLOT 


(BLXH-0. 1,0.1, 


2) 


CALL 


PLOT 


(0.1, 0.1, 2) 





TITLE 

CALL TITLE 
SCALE 

CALL METER( SCALE) 

READ THE DATA 
10 CONTINUE 

RE AD ( 51 , END =50 )MONTH , 0 AY , YE AR , HOUR , 

$ MIN , SECS , LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG 



JOB (0812,9999), ‘EZE(3UIEL' ,CLASS=C 
EXEC FRTVCLGP 
DD * 

PLOT 

FORTRAN VS 
AUGUSTO EZEQUIEL 



so 
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COMPUTES THE UTM COORDINATES 

CALL GPUTM( L ADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG , YPOS ,XPOS ) 
COMPLTTES THE TIME IN SECS 

TIME2=SECD*FL0AT( DAY )+SECH*FLOAT( HOUR )+60 . D05^FL0AT( MIN )+SECS 

COMPUTES THE PLOTTER COORDINATES 

XPLT=( XPOS-XLE FT )*SCALE5^100 . 

YPLT=( YPOS-YLEFT )*SCALE*100 . 

TESTE IF INSIDE AREA 

IF(XPLT.LT.O. ) GOTO 30 
IF( YPLT.LT.O. ) GOTO 30 
IF(XPLT.GT.BLXH) GOTO 30 
IF(YPLT.GT.BLYH) GOTO 30 

TESTS IF POSITIONS ARE AHAY MORE THEN 10 SECS IN TIME 

IF(TIME2-TIME1.GT.10.D0) IPEN=3 
TIME1=TIME2 

PLOTS THE POSITION 

CALL PLOT(XPLT,YPLT,IPEN) 

IPEN=2 

PLOTS THE TIME EVERY 05 MINUTES 

IF(MODtMIN,05KGT.O.OR.IFIX(SECS).GT.O) GO TO 10 
VALUE = FLOAT( HOUR )*100+FLOAT( MIN ) 

CALL NUMBER(XPLT+0.15,YPLT+0. 15,0. 25, VALUE ,0. ,-l) 

CALL SYMBOL(XPLT,YPLT, 0.15,3,0. ,-l) 

GOTO 10 

PEN UP WHILE THE POSITIONS ARE OUT OF THE SHEET 

30 IPEN=3 

GOTO 10 

END OF PLOT 

50 CALL PL0T( 0. ,0. ,+999) 

STOP 

END 

SUBROUTINE METER( SCALE) 

X0=24.0 

Y0=2.0 

CALL PL0T(X0,Y0,-3) 

DO 10 J=l,12 
XP= FL0AT(J)-1.0 
CALL PLOT(XP,0.0,3) 

CALL PL0T(XP,0.25,2) 

10 CONTINUE 

DO 20 J=l,9 
XP=FL0AT( J)*0.1 
CALL PL0T(XP,0.0,3) 

CALL PL0T(XP,0.2,2) 

20 CONTINUE 

CALL PL0T(0. ,0.0,3) 

CALL PLOTdl. ,0. ,2) 
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CALL PLOT(0. ,0.2,3 ) 

CALL PLOTdl. ,.2,2) 

C 

VALUE=1. /SCALE/100. 

CALL NU^BER( -0.25,0. 27, 0.25, VALUE, 0. 0,-1) 

CALL NUMBER( 0.95,0. 27,0. 25,0.0,0.0,-! ) 

VALUE=5.0/SCALE/100. 

CALL hA>1BER( 5. 775, 0.27, 0.25, VALUE, 0. 0,-1) 

VALUE=10 . O/SCALE/100 . 

CALL NUMBERC 10 . 7,0 . 27 ,0 . 25 ,VALUE ,0 . 0 , -1 ) 

CALL SYMBOL(5.5,-0.50,.25,6HMETERS,0.0,6) 

C 

CALL $YMBOL(5.5,24.0,2.0,62,0.0,-1) 

CALL SYMBOL! 5. 7, 2*+. 5,0. ^,85,0. 0,-1) 

XO=-XO 

YO=-YO 

CALL PL0T(X0,Y0,-3) 

C 

RETURN 

END 

C 

SUBROUTINE TITLE 
C 

X0=26.5 

Y0=43.05 

CALL PL0T(X0,Y0,-3) 

C 

CALL SYMBOL(00.0,3.5,0.25,29H NAVAL POSTGRADUATE SCHOOL ,0.,29) 
C 

CALL SYMBOL! 00. 0,3. 0,0. 25, 29HSEAFLOOR BENCHMARK EXPERIMENT ,0 ., 29 ) 



CALL SYMBOL! 00.0,2.5,0.25,29H 
CALL SYMBOL!00.0,2.0,0.25,29H 
CALL SYMBOL! 00.0,1.5,0. 25, 29H 
CALL SYMBOL! 00. 0,1. 0,0. 25, 29H 
CALL SYMB0L!00.0,1.0,0.25,29H 
CALL SYMBOL! 00. 0,1. 0,0. 25, 29H 
CALL SYMBOL!00.0,0.5,0.25,29H 
CALL SYMBOL! 00.0,0.5,0.25,29H 



PHASE II 
R/V POINT SUR 
16/17 AUGUST 1986 
GESAR SOLUTIONS 
TI<+100 

MINI RANGER FALCON 
ANTENNA POSITIONS 
TRANSDUCER POSITIONS 



,0.,29) 
,0. ,29) 
,0. ,29) 
,0. ,29) 

,0. ,29) 
,0. ,29) 

,0. ,29) 
,0. ,29) 



CALL SYMBOL!00.0,0.0,0.25,29H LCDR AUGUSTO EZEQUIEL ,0.,29) 

CALL PLOT!-0.2,-0.2,3) 

CALL PLOT! -0.2,3.95,2) 

CALL PLOT!7.5,3.95,2) 

CALL PLOT!7.5,-0.2,2) 

CALL PLOT! -0.2, -0.2, 2) 



XO=-XO 

YO=-YO 

CALL PLOT!XO,YO,-3) 

RETURN 

END 



C 



SUBROUTINE GPUTM! LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG , NORTH , EAST ) 
DOUBLE PRECISION A ,R ,N , AP ,BP ,CP ,DP ,EP ,S ,R1 ,ESQ , ESQP ,RM,RP ,KO 
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DOUBLE PRECISION R2 >R3 ,R4 ,R5,P ,P2 ,P3 ,P4,P5>P6 ,A6 ,B5 >SINSEC 
DOUBLE PRECISION PHI ,0LAM, NORTH , EAST ,B,PHIMIN, PI , LON 
REAL LASEG,LOSEG 

THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 123 00 00 W 

CM=-123.000 

PHI=DFLOAT ( LAOEG )+OFLOAT( LAMIN )/60 . D0+DBLE( LASEG 1/3600 . DO 
LON=DFLOAT ( LOOEG )+OFLOAT( LOMIN 1/60 . 00+DBLE( LOSEG 1/3600 . DO 
LON=-LON 

DLAM=( LON-CM 1*3600. ODO 

A=6378135.0D0 
R=298.26D0 

KO=0.9996D0 

B = A*(R-1.D0)/R 

N = (A-B1/(A+B1 

AP = A*( ( l.DO-N 1+5. DO/4. D0*(N**2-N**3 1+81. 00/64. 00*(N**4-N**51) 

BP = 3. 00/2. D0*A*( (N-N**21+7. DO/8. D0*(N**3-N**4)+55. 00/64. 00*N**51 
CP = 15.D0/16.00*A*(N**2-N**3+3.00/4.D0*(N**4-N**51 1 
DP = 35. 00/48. 00*A*( N**3-N**4+ll . DO/16 . D0*N**5 1 
EP = 315. DO/512. D0*A*(N**4-N**51 
PHIMIN = PHI*60. 00*2.908882086660-4 

PI=OARCOS( -l.DOl 

PHI=PHI/180.00*PI 

S = AP*PHIMIN-BP*0SIN( 2.00*PHI 1+CP*0SIN(4.00*PHI 1 
$ - 0P*0SIN( 6.00*PHI 1+EP*DSIN(8.00*PHI 1 

R1 = KO*S 

SINSEC = ( 1.00/3600.00 1/180. 00*PI 
SINSEC=0SIN(SINSEC1 

ESQ = ( A**2-B**21/A**2 
ESQP = ESQ/( 1. DO-ESQ 1 

RM = A*(l. DO-ESQ 1/(0SQRT(1.00-ESQ*DSIN( PHI 1**211**3 
RP = RM*( 1.00+ESQP*0COS( PHI 1**21 

R2 = RP*DSIN(PHI1*DCOS(PHI1*SINSEC**2/2.DO*KO*1.0D8 
R3 = SINSEC**4*RP*OSIN(PHI1*DCOS(PHI)**3/24.DO*(5.DO-OTAN(PHI1**2 
$ + 9.DO*ESQP*OCOS(PHI1**2+4.00*ESQP*ESQP*DCOS(PHI1**41*KO*1.D16 
R4 = RP*DC0S(PHI1*SINSEC*K0*1.D4 

R5 = SINSEC**3*RP*OCOS(PHI)**3/6.DO*(1.DO-OTAN(PHI1**2 
$ + ESQP*0C0S( PHI 1**2 )*K0*1. 012 

P = .0001D0*DLAM 
P2 = P**2 
P3 = P**3 
P4 = p**4 
P5 = P**5 
P6 = P**6 

A6 = P6*SINSEC**6*RP*0SIN(PHI 1*OCOS(PHI l**5/720.00 
$ * ( 61. 00-58. DO*DTAN( PHI 1**2+DTAN( PHI )**4 

$ + 270 . DO*ESQP*OCOS( PHI l**2-330 . DO*ESQP*DSIN( PHI 1**2 1*K0*1 . 024 

B5 = P5*SINSEC**5*RP*OCOS(PHI1**5/120.DO*(5.DO-18.DO*DTAN(PHI 1**2 
$ + DTAN(PHI1**4+14.00*ESQP*OCOS(PHI)**2 

$ - 58.D0*ESQP*0SIN( PHI 1**2 1*KO*1.D20 

N0RTH=R1+R2*P2+R3*P4+A6 
EAST=( R4*P+R5*P3+B5 1+500000 . DO 
RETURN 
END 

/* 

//GO.PLOTPARM DO * 
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8PL0T XMIN=0 
//GO.FT06F001 
//GO.FT51F001 

// 



. ,XMAX=S8. ,YMIN=0. ,YMAX=51. ,SCALE=1. ,UNITS=2.S40 SEND 
DD sysoi;t=» 

00 0SN=MSS . S0812 . GPS . POS .ANTENNA , OISP=SHR 
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APPENDIX E 

PROGRAM CVFICA. SOURCE LISTING 



//JOBCVFI JOB (0812,9999), 'EZEQUIEL' ,CLASS=C 
//*MAIN 0RG=NPGVM1 . 0812P ,LINES=( 99 ) ,CARDS=( 99 ) 

//*F0RMAT PR,DDNAME=GO.FT06F001, 

//*FORMS=SEPl 

// EXEC FORTVCLG 

//FORT.SYSIN DO * 

PROGRAM CVFICA 

C AUTHOR: AUGUSTO EZE(?UIEL 

C DATE: MARCH 17, 1987 

C 

C DESCRIPTION: 

C THIS PROGRAM CONVERTS THE FICA FILES INTO THE 

C INPUT FORMAT OF THE KALMN PROGRAM AND GIVES 

C GENERAL INFORMATION ABOUT THE DATA 
C 

DOUBLE PRECISION FPI ( 500 ) ,TP ,PB, HD ,RL10FF ,RL20FF 
DOUBLE PRECISION Q1 ,Q2 ,PI ,C ,VOS ,K1 ,K31 ,K32,TLDLL ,TLPLL 
DOUBLE PRECISION PDBPC ,PDBC ,CN1( ^ ) ,CN2( 4 ) ,CR1( ^ ) ,CR2( 4 ) 

DOUBLE PRECISION DOPK <+ ) ,DOP2( <+ ) ,BLC ,BPDC ,SGR1( 4 ) ,SGR2( 4 ) 

DOUBLE PRECISION SGDK 4 ) ,SGD2( <♦ ) ,K2,K4,TT ,M,AG,ES(5 
DOUBLE PRECISION TOCS,AS( 3 ) , ADC ,CRSS ,DNS,MOS,CUCS,ES,CUSS 
DOUBLE PRECISION SQAS ,TOES ,CICS ,OMEGS,CISS,IOS,CRCS ,WS 
DOUBLE PRECISION DLAT ,DLOG,DHT 
DOUBLE PRECISION OMEDS ,ADE ,IDOTS,X,Y,Z,ROBL 
REAL LASEC,LOSEC,SEC 

INTEGER INTG( 500), BLOCK, TYPE, NCI, Nil, NFI,LAMIN,LOMIN,LADEG,LODEG 
INTEGER TRACER, MQVEL( 4 ),ISAT( 4), INDEX, SAT, IWKNOS,IEDATS,NSAT 
INTEGER 10, STDAY, MONTH, YEAR, HOUR, MIN, DAY 
INTEGER STAT(4) 

CHARACTER*8 CHK500 ),TITLE( 10 ) 

CHARACTER*5 BLK 
LOGICAL SAVE, FIRST 
C 

C SETS DEFAULT METEO DATA 

C 

DATA TP,PB,HD /15 . , 980 . , 75 ./ 

C 

C SETS THE PARAMETERS OF HGS72 

C 

DATA AG, ROBL /6378 . 1350D0 , 298 . 26D0/ 

C 

C FLAG TO SAVE DATA ONLY HHEN THERE IS NAVIGATION DATA 
C FOR ALL TRACERS 

C 

SAVES. FALSE. 

C 

C FLAG TO SAVE INFORMATION DATA AS TYPE 1 ONLY ON FIRST 
C CASSETTE 

C 

FIRSTS. TRUE. 

C 

C ECCENTRICITY SQUARED OF REFERENCE ELLIPSOID 
C 

ESQ=( 2 . D0-( 1 . DO/ROBL ) )/ROBL 
C 
C 

c 

Q1=154.D0 

Q2=120.D0 

C 

C COMPUTS THE VALUE OF PI IN DOUBLE PRECISION 
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C 



PI=DARCOS( -1.00) 



SPEED OF LIGHT (KM/SEC) 

0=299792.^5800 
NOMINAL SATE LITE FREQUENCY 
VOS=10.2306 

CALCULATION OF CONSTANTS TO BE USED IN THE SGIMA COMPUTATIONS 
K1=(CA'0S)**2 

K31 = ( C/( 2 . DO^PIWOS^Ql ) )**2 
K32= ( C/( 2 . D0*PI W0S*Q2 ) )**2 



RL10FF=-6000.DO 

RL20FF=7600.DO 



M=4.D0 

TLDLL=.7D0 

TLPLL=.7D0 

PDBPC=100.D0 

PDBC=100.00 

SETS THE FLAGS FOR NOT HAVING THE EPHEMIS DATA 

DO 3 1=1,4 
3 STAT(I)=0 

PRINTS OUT THE IDENTIFICATION OF THE PROGRAM 



WRITE! 6,2) 

FORMAT! IHl ,///////// 



,/> • 

,/, ' 

>/> • 

,/,• 

>/> ' 

,/////////) 



OUTPUT OF THE PROGRAM CVFICA THAT CONVERTS ' 
** THE FICA FILES TO A FORMAT READABLE BY THE ‘ 
PROGRAM KALMN2 

' 



READS THE TITLE FROM THE FILE AND CONTROL INFORMATION 

READ!31,5,END=2000)TITLE 
WRITE! 6,5) TITLE 
READ! 31,5,END = 2000 )TITLE 
WRITE! 6,5, END=2000) TITLE 
5 FORMAT! 10A8) 

READ! 31,*,END = 2000 ) MONTH ,STDAY,YEAR 
WRITE! 6,6) MONTH, STDAY, YEAR 
6 FORMAT!///, 

^ ■ DATE OF STARTING GPS WEEK ! MONTH DAY YEAR) ',13,13,15) 

INPUTS THE DATA FROM THE COMMAND FILE 

READ! 31,*,END=2000) NSAT, K2, K4 
READ!31,*,END=2000) X,Y,Z 
WRITE! 6,1) NSAT,K2,K4,X,Y,Z 
1 FORMAT ! ' NUMBER OF TRACERS ',120,/, 

n ' RANGE SIGMA BIAS FACTOR ' ,E20.8,/, 
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« • DOPPLER SIGMA BIAS FACTOR ',E20.8,/, 

?? ' ESTIMATE RECEIVER POSITION',/, 

^ ' X COORD' ,F20.6,/, 

^ ‘ Y COORD' ,F20.6,/, 

^ ' Z COORD' ,F20.6,/) 

C 




c 

C READS THE CONTROL BLOCK INFORMATION 
C 



10 READ( 50,20,END=2000 ) BLK , BLOCK ,NFI ,NII ,NCI 
20 FORMAT! AS, 415) 

C 

C READS THE FLOATING PART DATA IF ANY 
C 

IF (NFI.EQ.O) GOTO 35 

READ! 50,30, END = 2000 ) ! FPI ! I ) ,1=1 ,NFI ) 

30 FORMAT! 4D20. 16) 

C 

C READS THE INTEGER PART IF ANY 

C 

35 CONTINUE 

IF !NII.EQ.O) GOTO 45 

READ! 50,40 ,END=2000 ) ! INTG! I ) ,1=1 ,NII ) 
40 FORMAT! 6112) 

C 

C READS CHARACTER PART DATA IF ANY 

C 

45 CONTINUE 

IF !^K:I.EQ.O) GOTO 55 
READ!50,50,END=2000) !CHI!I ), 1 = 1, NCI) 

50 FORMAT! 10! AS)) 

55 CONTINUE 



IF !BL0CK.NE.6) GOTO 100 




C 

C USER EPOCH TIME OF PSEUDO RANGE 
C 

TT=FPI!3) 

C 

C COMPUTES THE DATE 

C 

CALL SUB! TT, DAY, HOUR, MIN, SEC) 

DAY=DAY+STDAY 

C 

C IF NO EPHMERIS DATA AVOIDS COMPLTTATIONS 
C 

IF! .NOT. SAVE) GOTO 10 
C 

C LI CARRIER SIGNAL TO NOISE 

C 

CN1(1)=FPI!4) 

CN1(2)=FPI!S) 

CN1(3)=FPI!6) 

CN1!4)=FPI! 7) 

C 

C L2 CARRIER SIGNAL TO NOISE 

C 

CN2(1)=FPI!8) 

CN2!2)=FPI!9) 

CN2!3) = FPI! 10) 

CN2(4)=FPI!11) 

C 

C LI PSEUDO RANGE !KM IN FICA FILES CONVERTED TO SECONDS FOR KALMN2 ) 
C 

CR1!1)=FPI!12)/C 
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CRK 2) = FPI( 13)/C 
CR1(3)=FPI(14)/C 
CR1(^)=FPI(1S)/C 

L2 PSEUDO RANGE (KM IN FICA FILES CONVERTED TO SECONDS FOR KALMN2 ) 

CR2(1)=FPI(16)/C 
CR2( 2) = FPI( 17)/C 
CR2(S) = FPI( 18 VC 
CR2(^)=FPI( 191/C 

LI CARRIER DOPLER PHASE 

DOPK 1) = FPI( 20) 

D0P1(2)=FPI(21) 

D0P1(3)=FPI(22) 

D0Pl(^f) = FPI( 23) 

L2 CARRIER DOPLER PHASE 

D0P2(1) = FPI(2<^) 

D0P2(2)=FPI(25) 

D0P2(3)=FPI(26) 

D0P2(4) = FPI( 27) 

SV PRN OF EACH TRACER 

ISAT(1)=INTG(1) 

ISAT(2)=INTG(2) 

ISAT(3)=INTG(3) 

ISAT(4)=INTG(^) 

L1,L2 QUALITY FACTOR ( TRACER , FREQUENCY ) 

MQVEL( 1)=INTG( 9)-INTG( 13) 

MQVEL( 2)=INTG( 10)-INTG( 14) 

MQVEL( 3 )=INTG( 11 )-lNTG( 15 ) 

MQVEL( 4 ) = INTG( 12 )-INTG( 16 ) 

IF ANY ERRORS FLAGS THE DATA OF THE SATELLITE 

DO 120 INDEX=1,NSAT 

ERRORS IN SIGNAL 

IF( MQVEL( INDEX ).EQ.O) GOTO 114 
MQVEL( INDEX )=15 
I1=8+INDEX 
I2=12+INDEX 

WRITE (6,110) MONTH , DAY , YE AR , HOUR , MIN , SEC , INDEX , INTG( 1 1 ) ,INTG( 1 2 ) 
110 FORMAT(/,’ WARNING AT 

« 2(12, IX), 14, IX, 2(12,': ' ),F5.2,/, 

» ■ DATA FROM TRACER ’ , 

U 12, • GESAR BAD STATUS ‘,15, IX, IS) 

MQVEL( INDEX) =15 

IF THE DATA IS ALREADY BAD NO NEED TO TEST THE SATELLITE 
GOTO 120 

NON EXISTENT SATELITE DATA FOR THAT TRACER 
114 CONTINUE 

IF( STAT( INDEX ).EQ.ISAT( INDEX)) GO TO 120 



WRITE (6,115) MONTH , DAY ,YEAR , HOUR , MIN , SEC , INDEX, 
n ISAT( INDEX ) ,STAT( INDEX ) 

115 FORMAT (/, • ^ WARNING AT ’ , 
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» 2X,I2,1X,I2,1X,I^,2X,I2, • : • ,12 , • : • , F5 . 2 ,/ , 

^ ‘ TRACKING TRACER ',I2,‘ SV ‘,12, 

EXISTING NAV DATA FOR SV M2,‘ FLAGEO BAD QUALITY DATA') 
MQVEU INDEX ) = 15 
12D CONTINUE 
C 

C COMPUTS THE SIGMAS 

C 

DO 15D INDEX=1,NSAT 
IF( MQVEL(INDEX).NE.D) GOTO 15D 
C 

C CALCULATE LI AND L2 SIGMA RANGE 

C 

BLC=TLDLL 

BPDC=PDBPC 

SGRK INDEX ) = ( (BLC)^)/( ID. CNK INDEX )/lD . ) ) ) 

» * ( .55^BPDC/aD.*^^{CNl( INDEX )/lD. )) )i^l+K2 

SGR2( INDEX ) = ( (BLC*M)/( 1D.**( CN2( INDEX )/lD. ) ) ) 

« * ( .5*BPDC/( ID. **(CN2( INDEX )/lD. )))*K1+K2 

C 

C CALCULATE LI AND L2 SIGMA DOPPLER DATA 
C 

BLC=TLPLL 

BPDC=PDBC 

SGDK INDEX ) = ( ( BLC*M )/( ID CNK INDEX )/lD . ) ) ) 

^ * (l.+( .5*BPDC/(10.i<^(CNl(INDEX)/lD. ))))*K314K<+ 

SGD2( INDEX ) = ( ( BLC^^M )/( ID . CN2( INDEX )/lD . ) ) ) 

^ * ( l.+( .5*6PDC/( 1D.**(CN2( INDEX)/1D. ) ) ) )*K32+K<^ 

15D CONTINUE 
C 

C OUTPUT THE DATA TO THE FILE IF THERE IS NAV DATA 
C 

TYPE=3 

WRITE! 7) TYPE 

WRITE! 7) TT, ISAT , CRl, CR2, DOPl, DOP2, 

^ SGRl, SGR2, SGDl, SGD2, MQVEL 

C 

C GO TO READ ANOTHER BLOCK 
C 



GOTO ID 



IDD CONTINUE 

IF (BLOCK. NE.1D9) GOTO 2DD 



C^* NAVIGATION DATA ! AS TRANSMITED) 






C 



C 

C 

c 

c 

2DD 



C 



TRACER=INTG! 1) 

SAT=INTG( 2) 

GO TO READ ANOTHER BLOCK 
GOTO ID 
CONTINUE 

IF !BLOCK.NE.9) GOTO 3DD 



NAVIGATION DATA ! DECODED SUBFRAMES 1 TO 3 ) 






C 

C DAY OF WEEK 

C 

IWKN0S=DINT!FPI!6)) 

C** IWKNOS=IAND( IWKNOS , 1777B ) 

C 

C SV HEALTH 

C 

IEDATS=DINT!FPI!9)) 



C 
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CLOCK EPOCH (GPS SECONDS OF WEEK) 

T0CS=FPI(13) 

CLOCK BIAS (SEC) 

AS(1)=FPI(16) 

CLOCK DRIFT (SEC/SEC) 

AS( 2 ) = FPI( 15) 

CLOCK DRIFT RATE ( SEC/SEC»*2 ) 

AS(5)=FPI(1^) 

AGE OF DATA (CLOCK) (SEC) 

ADC=FPI(10) 

RADIAL SINE CORRECTION (DIVIDED BY 1000. TO GET KM) 
CRSS=FPI( 27)/1000. 

CORRECTION TO MEAN MOTION (RADIANS/SECONDS) 
0NS=FPI(28) 

MEAN ANOMALY AT EPOCH (RADIANS) 

MOS=FPI( 29) 

IN TRACK COSINE AMPLITUDE (RADIANS) 

CUCS=FPI(50) 

ECCENTRICITY 
ES=FPI(31) 

IN TRACK SINE AMPLITUDE (RADIANS) 

CUSS=FPI(32) 

SQUARE ROOT OF SEMI -MAJOR AXIS CONVERTED TO SQR OF KM 

SQAS=FPI(33 
SQAS=SQAS/1000.D0 
SQAS=DSQRT(S(3AS) 

TIME OF EPOCH (GPS SECONDS OF WEEK) 

T0ES=FPI(34) 

I^K:LINATION COSINE CORRECTION (RADIANS) 

CICSaFPI(<^6) 

RIGHT ASCENSION NODE (RADIANS) 

0MEGS=FPI(47) 

INCLINATION SINE CORRECTION (RADIANS) 

CISS=FPI(^8) 

INCLINATION ( RADIANS ) 

I0S=FPI(<^9) 

C 
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c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 



RADIAL COSINE ADJUSTMENT (DIVIDED BY 1000 TO GET KM ) 
CRCS=FPI(SO )/1000. 

ARGUMENT OF PERIGEE (RADIANS) 

WS=FPI(51) 

RIGHT ASCENSION OF ASCENDING NODE (TIME DERIVATIVE) 

( RADIANS/SEC ) 

OMEDS=FPI(52) 

AGE OF DATA (SEC) 

ADE=FPI( 26) 

INCLINATION TIME DERIVATIVE (RADIANS/SECONDS) 
ID0TS=FPI(54) 

SETS THE FLAG FOR EXISTANCE OF NAV DATA FOR THIS SATELITTE 



ST AT( TRACER )=SAT 
C 

C TESTS IF THERE IS ENOUGH NAV DATA 
C 

SAVE=.TRUE. 

DO 210 INDEX=1,NSAT 
I F ( ST AT ( IND EX ) . EQ . 0 ) SAVE = . FALSE . 

210 CONTINUE 
C 

C OUTPUTS THE DATA 
C 

TYPE=2 

WRITE(7) TYPE 

WRITE ( 7 ) TRACER, SAT, AS( 1 ) ,AS( 2 ) ,AS( 3 ) ,CICS,CISS 
n ,CRCS,CRSS, cues, cuss, DNS, ES,IDOTS,IOS,MOS 

n , OME OS , OME GS , SQAS , AO E , TOCS , TOES , I WKNOS , WS 

n ,IEDATS,AOC 

c 

C INFORMS THAT NAVIGATION DATA WAS RECEIVED FOR THIS SATELLITE 

C 

C 

WRITE (6,220) TRACER , SAT , MONTH , DAY , YEAR , HOUR , MIN , SEC 
220 FORMAT (///,' ********^HHHHt*****************)^*********** • ,/ 

n ,' ** NAVIGATIONAL **',/, 

n • ** DATA FROM TRACER ',12,' SV ',13,' RECEIVED **',/, 

n ' • ,7X,I2,1X,I2,1X,I<^,2X,I2, • : • ,12, • : ' ,F5.2,7X, 

c 

C PRINTS OUT THE WARNING WHEN THE HEALTH IS NOT GOOD 
C 



IF( FPK9) .NE. 0.0 ) WRITE(6,201) TRACER , SAT ,FPI ( 9 ) 
TRACER ',13,' SV ',13,' HEALTH STATUS ' ,F15.6, 

ff p ^ A A ^ A A A! A A A W K A X ^ K A A A AAA | 

C 

C GO TO READ ANOTHER BLOCK 
C 

GOTO 10 
C 

300 CONTINUE 

IF ( BLOCK. NE. 101) GOTO 400 
C 



C** BLOCK 1 PLUS INPUT DATA 









c 

C PRESSURE 
C 
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PB=FPI(1) 



C 

C TEMPERATURE 
C 

TP=FPI( 2) 



C 

C HUMIDITY 

C 

HD=FPI(3) 



C 

C TESTE FOR BAD WEATHER DATA, CHANGE TO STANDARD VALUES BAD DATA 
C 



IF 


( ( 


PB 


.LT. 


800. 


) 


.OR. 


( PB 


.GT. 


1200. 


) ) 


PB = 


980 


IF 


( ( 


TP 


.LT. 


-99. 


) 


.OR. 


( TP 


.GT. 


+ 99. 


)) 


TP = 


15 


IF 


( ( 


HD 


.LT. 


1. 


) 


.OR. 


( HD 


.GT. 


100. 


) ) 


HD = 


75 



C 

C OUTPUTS THE DATA FOR THE FIRST TIME 
C 

IF( .NOT. FIRST) GOTO 10 
C 

TYPE=1 

WRITE(7) TYPE 

WRITE( 7) RL10FF,RL20FF,PB,TP,HD,X,Y,Z 
C 

C SETS THE FLAG 

C 

FIRSTS. FALSE. 

C 

C GO TO READ ANOTHER BLOCK 

C 

GOTO 10 
C 

<400 CONTINUE 

IF (BLOCK. NE. 5) GOTO 460 
C 



C** GESAR SOLUTION 
C 

C SOLUTION USER EPOCH TIME 
C 

TT=FPI(2) 

C 

C IF NO EPHMERIS DATA DOES NOT KEEP THE DATA 
C 

IF( .NOT. SAVE) GOTO 10 
C 

C COORDINATES OF RECEIVER 

C 

X=FPI(3) 

Y=FPI(4) 

Z=FPI(5) 

C 

C COMPUTES THE DATA 

C 

CALL SUB(TT, DAY, HOUR, MIN, SEC) 

DAY=DAY+STDAY 

C 

C ONLY PRINTS POSITIONS IF UNCOMENTED 

C 

C WRITE(6,450) MONTH , DAY , YEAR , HOUR ,MIN,SEC ,X,Y ,Z 

450 FORMAT (/, ' GESAR SOLUTION AT (MM DAY YYYY HH :MM:SS. SSS ) ' , 

^ 2X,I2,1X,I2,1X,I4,2X,I2, ' ; ' ,12, ' : * ,F5.2, 

»/,’ PSN COORD X, Y, Z (KM) ',F14.3,' ',F14.3,' ',F14.3) 

C 

C CONVERTS TO LATITUDE / LATITUDE / HT 

C 

CALL XYLLH(X,Y,Z,AG,ESQ,DLAT,DLOG,DHT,PI ) 

C 

C PRINTS THE RESULTS 
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c 

C WRITE( 6,455) DLAT ,DLOG ,DHT 

455 FORMAT ( • LAT, LON, HT (DEG, DEG, KM) ’ ,2F17. 12 ,2X,F12. 7 ) 

C 

C CONVERTS LATITUDE TO DEGREES MINLTTES AND SECONDS 
C 

CALL DMS( DLAT,LADEG,LAMIN,LASEC) 

C 

C CONVERTS LONGITUDE TO DEGREES MINITTES AND SECONDS 
C IN A RANGE 0 TO 180 

IF (DLOG.GT,180.DO) DL0G=360 . DO-DLOG 
CALL DMS(DLOG,LODEG,LOMIN,LOSEC) 

C 

C WRITE (6,456) LADEG , LAMIN , LASEC , LODEG , LOMIN , LOSEC 

456 FORMAT ( • LAT, LON (DEG MIN SEC, DEG MIN SEG ) • ,2( 2I4,1X,F7.4 ) ,/ ) 
C 

C OUTPUTS THE DATA TO A FILE 
C 

WRITE ( 8 ) MONTH , DAY , YEAR , HOUR ,MIN ,SEC , 

» LADEG , LAMIN , LASEC , LODEG , LOMIN , LOSEC 

C 
C 

C GO TO READ ANOTHER BLOCK 
C 

GOTO 10 
C 

460 CONTINUE 

IF ( BLOCK. NE. 13) GOTO 470 
C 

C^^* TAPE HEADER/TRAILOR ** 



C 



461 



462 

463 

464 

465 

466 



WRITE( 6,461 ) 

FORMAT ( / , 5X , ' x , 

n 5X,'x* TAPE HEADER/TRAILOR **’,/, 

tr OA > ^ ^ y/ J 

IF( INTG(l).EQ.l) WRITE (6,462) 

IF( INTG( 1).E(3.2) WRITE (6,463) 

IF( INTG( D.EQ.4) WRITE (6,464) 

IF( INTG(1).EQ.8) WRITE (6,465) 

FORMAT( • BEGINING OF DATA SET',/) 

FORMAT( • END OF DATA SET',/) 

FORMAT ( ' BEGINING OF CASSETTE',/) 

FORMAT ( ' END OF CASSETTE',/) 

WRITE( 6,466) INTG( 2 ) ,INTG( 3 ) 

FORMAT ( ' CASSETTE ',15,' IN DATA SET',/, 

U ' CASSETTE SEQUENCE NUMBER ',15) 



C 

C GO TO READ ANOTHER BLOCK 
C 

GOTO 10 
C 

470 CONTINUE 

IF (BLOCK.NE.il) GOTO 480 
C 

C IF NO EPHMERIS DATA DOES NOT KEEP THE DATA SO NO NEED TO 
C INFORM THE USER 

C 

IF( .NOT. SAVE ) GOTO 10 



C 



C5f* RECEIVER ERROR BLOCK 



C 

C 

WRITE(6,471) MONTH, DAY, YEAR, HOUR, MIN, SEC 

471 FO^RMAK/, ^ ^ ERROR BLOCK RECEIVED AFTER : **',/> 

;t ' ,9X,2(I2,1X),I4,2X,2(I2,': ' ),F5.2,8X,'**' ,/, 
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c 

HRITE (6,^72)aNTG(I),I = l,NII) 



FORMAT! ' 


ERROR LOG MESSAGE TI ICD FORMAT » X' 


,18,/, 


« ' 




2* 


,18,/, 


' 


n 


3' 


,18 ,/ , 






4' 


,18,/, 


n 


n 


5' 


,18,/, 


» 


u 


6' 


,18,/, 


n 


» 


7' 


,18,/, 






8' 


,18,/, 


n ' 


ERROR LOG OVERFLOW COUNT 


1 


,18,/) 



C 

C GO TO READ ANOTHER BLOCK 
C 

GOTO 10 
C 

480 CONTINUE 

IF (BL0CK.NE.8) GOTO 500 
C 

C IF NO EPHMERIS DATA DOES NOT KEEP THE DATA SO NO NEED TO 
C INFORM THE USER 

C 

IF( .NOT. SAVE ) GOTO 10 



C 




C 

C 

C 



C 

C 

481 



C 



482 



C 

C 

C 



490 



495 

C 

C 

C 



COMPUTES THE DATE 

CALL SUB(FPI(2),OAY,HOUR,MIN,SEC) 
DAY=DAY+STOAY 



WRITE! 6 ,481 ) MONTH , DAY , YEAR , HOUR, MIN ,SEC 

n ' TRACKING CONFIGURATION AT 

« ' ,9X,2(I2,1X),I4,2X,2(I2, • : ' ),F5.2,8X,'**' ,/, 



WRITE (6,482) FPI ( 1 ) , ( FPI ( I ) , I =3 ,NFI ) , ( INTG( I) , 1 = 1 ,NII ) 

FORMAT!' PSEUDORANGE FTF OF VALIDITY SEC ',E22.14,/, 

3 ’ PREDETECTION BANDWITHS HZ FREQUENCY ,TRACER ’,/ , 

» 2(4E20. 14, /),//, 

» • TRACKING LOOP BANDWITHS HZ FREQUENCY , TRACER ',/ , 

» 4( 4E20.14,/),/, 

n • LOOP ROUND TRIP CALIBRACION DELAYS IN SECONDS * ,/, 

INCLUDING PRE-ANTENNA* ,/,' LI *,E20.14,' L2 *,E20.14,//, 

^ • LOOP ROUND TRIP CALIBRACION DELAYS IN SECONDS ' ,/, 

' (INTERNAL TO RFM)',/,' LI ',E20.14,' L2 ',E20.14,//, 

■ SV PRN ID OF EACH TRACER ' ,/,4I15,//, 
n ' ANTENA INDICATOR/TRACER' ,/, 41 15,//, 
n * CODE INDICATOR FLAG/TRACER ',/ ,4115 ,// ) 

TESTE IF THE CONSTELATION MATCHES THE EXISTING EPHMERIS 
DATA 

DO 495 INDEX=1,NSAT 

IFISTAT! INDEX). EQ.INTG! INDEX)) GO TO 495 
WRITE! 6 ,490 ) INDEX, ST AT! INDEX ) ,INTG( INDEX ) 

FORMAT!////,' ******* WARNING iHHHf**',//, 

^ ' ACTUAL NAVIGATIONAL DATA FOR TRACER *,I2,' IS FROM SV ',12,/, 
n ' TRACKING CONFIGURATION SHOWS SV ',12,//) 

CONTINUE 

GO TO READ ANOTHER BLOCK 
GOTO 10 



C 

500 CONTINUE 
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BLCXKS THAT THIS PROGRAM IS NOT PREPERED TO DECODE ** 

C 



WRITE( 6,501) 

» /,’ ** FOUND A BLOCK THAT THIS PROGRAM IS NOT ** ' , 

^ ABLE TO DECODE . BLOCK DATA IS DUMPED NEXT , 

U / , • ^Ht***********i^**************^**********^*********** ' ) 

WRITE ( 6,520 ) BLK , BLOCK ,NFI ,NII ,NCI 
IF (NFI.NE.O) WRITE(6,530) ( FPK I ) , 1 = 1 ,NFI ) 

IF (NII.NE.O) WRITE(6,540) ( INTG( I ) ,1=1 ,NII ) 

IF (NCI.NE.O) WRITE(6,550) ( CHK I) ,1=1 ,NCI ) 

520 FORMAT( 1X,A5,«5) 

530 FORMAT! lX,‘^E20.1<h) 

5^0 FORMAT! IX, 61 12) 

550 FORMAT! IX, 10 !A8)) 

C 

C GO TO READ ANOTHER BLOCK 
C 

GOTO 10 
C 



C 

C END OF PROGRAM 
C 

2000 CONTINUE 
STOP 
END 
C 



C THIS SUBROUTINE WAS TAKEN FROM THE KALMN PROGRAM 

^ N^WWVfWWVfWWWWV^V^WV^V^WWV^^w y^N^WWV/V^WWV^ VrW WWV^WWWWV^WS^W. 

^ ^ ^ ^X ^X ^X ^X 7x ^X ^x ^X^xTx^X^X^X ^X ^X^X^X ^x 

c 

SUBROUTINE XYLLH 

G ! X , Y , Z , AG , ESQ , 

R DLAT , DLOG , DHT,PI) 

C THIS SUBROUTINE CONVERTS POSITION FROM THE EARTH-CENTERED CARTESIAN 
C COORDINATE SYSTEM TO THE GEODETIC COORDINATE SYSTEM. 

C X,Y,Z=ESTIMATE OF RECEIVER POSITION IN KM 
C AG=SEMIMAJOR AXIS OF REF. ELLIPSOID! RADIUS OF EARTH) IN KM 
C ESQ=ECCENTRICITY SQUARED OF REF. ELLIPSOID=! 2 . -1 ./OBL )/0BL 
C DLAT=GEODETIC LATITUDE OF RECEIVER POSITION IN DEGREES 
C DLOG=GEODETIC LONGITUDE OF RECEIVER POSITION IN DEGREES 
C DHT=GEODETIC HEIGHT OF RECEIVER POSITION IN METERS 
C PI VALUE OF PI COMPUTED IN DOUBLE PRECISION IN THE MAIN PROGRAM 
C SINCE LAT=LAT!X,Y,Z,LAT) AND HT=HT! X ,Y , LAT ) , AN ITERATIVE PROCEDURE 
C IS NECESSARY TO DETERMINE THE VALUES OF LATITUDE AND HT. 

C THE WORKING EQUATIONS ARE AS FOLLOWS: 

C HT=!R/C0S!LAT))-AG1 WHERE AG1=AG/! SQRT! 1-ESQ*SIN! LAT )**2 ) ) 

C TAN!LOG)=Y/X 

C TAN! LAT ) = ! ! Z+! ESQ^SIN! LAT )*AG1 )/R ) 

C 

DOUBLE PRECISION LT , LTl , LG , R , EP ,AG1 

DOUBLE PRECISION X , Y , Z > AG , ESQ , DLAT , DLOG , DHT , PI 
C 

C FIND LONGITUDE 

LG = DATAN2!Y,X) 

IF (LG .LT. O.DO) LG = LG ♦ !2. * PI) 

R=DSQRT!X**2+Y**2) 



C 

C A FIRST GUESS FOR LATITUDE WOULD BE ATAN!Z/R) 
LT1=DATAN!Z/R) 

C 

C SOLVE FOR LATITUDE BY ITERATIONS 

DO 10 1=1,5 

AG1=AG/! DSQRT! 1 . DO-ESQ*DSIN! LTl )**2 ) ) 
LT=DATAN! ! Z+! ESQ^DSIN! LTl )*AG1 ) )/R ) 
EP=DABS!LT-LT1) 

IF!EP.LE. .OOOOOOOOIDO ) GOTO 5 
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LTl = LT 
10 CONTINUE 

5 AGl=AG/( DSQRT( 1D0-ESG*DSIN( LTl )**2 ) ) 

C 

C FIND HEIGHT 

DHT=((R/DC0S(LT))-AG1) 

C CONVERT DLAT AND DLOG FROM RADIANS TO DEGREES 
DLAT=LT*180. DO/PI 
DLOG=LG*130. DO/PI 
RETURN 
END 
C 

C THIS SUBROUTINE CONVERTS DECIMAL DEGREES IN DEGREES MINUTES SECONDS 



c 

c 

c 



SUBROUTINE DMS( DEC ,LDEG ,MIN ,SEC ) 
DOUBLE PRECISION DEC ,XNUM ,XMIN 



XNUM=DABS( DEC ) 

LDEG=DINT(XNUM) 

XMIN=( XNUM-DFLOAT ( LDEG ) )*60 . DO 
MIN=DINT(XMIN) 

XNUM=( XMIN-DFLOAT( MIN ) )*60 . DO 
SEC=SNGL(XNUM) 

IF ( LDEG, GE. 360) LDEG=LDEG-360 

RETURN 

END 



C THIS SUBROUTINE CONVERTS TIME TAG IN DAYS HOURS MINUTES SECONDS 






C 



SUBROUTINE SUB ( TT , D AY , HOUR , MIN ,SEC ) 

INTEGER DAY, HOUR, MIN 
DOUBLE PRECISION TIME ,TT 
REAL SEC 
TIME=TT 

DAY = DINK TIME / ( 2^ . D0*3600 . DO ) ) 

TIME = TIME - DFLOAKDAY) * 2^. DO * 3600, DO 
HOUR = DINK TIME / 3600. DO) 

TIME = TIME - DFLOAKHOUR) * 3600, DO 
MIN = DINK TIME / 60, DO) 

TIME = TIME - DFLOAKMIN) * 60. DO 
SEC = SNGL(TIME) 

RETURN 

END 



/* 




//GO.FT06F001 


DD 


//GO.FT07F001 


DD 


//GO.FT08F001 


DD 


//GO.FT31F001 


DD 


//GO.FT50F001 


DD 


// 


DD 


// 


DD 


// 


DO 


// 


DD 


// 


DD 


// 


DD 


// 


DD 


// 


DD 


// 





SYSOUT=* 

DSN=MSS , S081 2 . CVFICA , OUTPUT , DISP=SHR 
DSN=MSS , S081 2 . GESAR , POS , DISP =SHR 
DSN=MSS,S0812,CVDATA,DISP=SHR 
DSN=MSS . S0812 , GPS , TAPES . DATA! DL285 ) , DISP=SHR 
DSN=MSS , S081 2 . GPS , TAPES , DATA( DL28^f ) , DISP=SHR 
DSN=MSS , S0812 , GPS , TAPES . DATA! DL286 ) , DISP=SHR 
DSN=MSS , S0812 , GPS , TAPES , DATA ( DL287 ) ,DISP=SHR 
DSN=MSS.S0812, GPS, TAPES. DATA! DL288 ) ,DISP=SHR 
DSN=MSS. S0812 . GPS . TAPES . DATA( DL289 ) , DISP=SHR 
DSN=MSS,S0812. GPS. TAPES. DATA! DL290 ) ,DISP=SHR 
DSN-MSS. S0812. GPS. TAPES. DATA(DL291), DISP =SHR 
DSN=MSS , S0812 . GPS . TAPES . DATA( DL292 ) ,DISP=SHR 
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APPENDIX F 

PR0GR.4M TRANSDUC GPS. SOURCE LISTING 



//JOBTRANS JOB ( D812 , 9999 ) , ‘ E2EQUIEL ' ,CLASS=C 
//*MAIN ORG=NPGVM1.0812P,LINES=( 99),CARDS=( 99 J 

//^FORMAT PR,DDNAME=GO.FTD6FD01, 

//*FORMS=SEPl 

// EXEC FORTVCLG 

//FORT.SYSIN OD * 

C 

C PROGRAM TRANSDUC GPS 
C 

C AUTHOR: AUGUSTO E2EQUIEL 

C DATE: MARCH 3D, 1987 

C 

C DESCRIPTION: 

C 

C THIS PROGRAMS TAKES EACH POSITION OF THE ANTENNA, 

C THE COMPUTED COURSE ( FILTERED OR NOT ) ,THE PITCH 

C AND ROLL DATA, AND COMPUTES THE POSITIONS OF THE 

C TRANSDUCER, APLYING THE SEVEN PARAMETER TRANSFORMATION 

C TO THE OFFSETS OF THE TRANSDUCER IN RELATION TO THE ANTENNA 

C CONSIDERED THE CENTER OF COORDINATE SYSTEM AND THEN ADDING 

C THE CORRECTIONS TO THE COORDINATES. 

C 

C THE PROGRAM RUNS IN MVS 

C 

C I/O SPECIFICATIONS: SEE END OF THIS JOB. 

C 

C THIS PROGRAM IS RESTRICT AS IS TO THE SEAFLOOR BENCHMARK 
C EXPERIMENT 

C 

C THE PROGRAM HILL RUN HITH ANY AMOUNT OF DATA ONLY LIMITED TO DISK 
C SPACE 

C 

C ANY BLANK LINES HILL TERMINATE THE PROGRAM IN ERROR 
C 
C 
C 

DOUBLE PRECISION XPOS,YPOS, TIMEl ,TIME2 ,TIME3 ,SECD ,SECH ,RATE 

DOUBLE PRECISION TIME^, TIMES 

DOUBLE PRECISION TT ,DLAT,DLOG,DHT ,X,Y ,2 

REAL COURSE,PITCH( 2),ROLU 2),SEC1( 2),SEC2,PIT,ROL,SEC3(2) 

REAL DX,DY,02,OFFX,OFFY,OFF2,HEAD( 2),LASEG,LOSEG 
INTEGER MONTHK 2 ) ,DAY1( 2 ) ,YEAR1( 2 ) ,HOURl( 2 ) ,MIN1( 2 ) 

INTEGER MONTH2 ,DAY2 ,YEAR2 ,HOUR2 ,MIN2 
INTEGER MONTH3(2),DAY3( 2),YEAR3( 2),HOUR3(2),MIN3(2) 

INTEGER LADEG,LAMIN,LODEG,LOMIN, STDAY , NSAT 
CHARACTER'S TITLE! ID) 

C 

C INITIALIZATION OF CONSTANTS 

C 

SECH*3600.DD 

SEC0=24.DD*SECH 

C 

0FFX=-4.268 

OFFY=-2.9SS 

OFF2=-9.502 

C 

C READS THE TITLE FROM THE FILE 

C 

READ(3D,S,EN0=10D)TITLE 
HRITE(6,5) TITLE 
REA0(3D,5,EMD=10D)TITLE 
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HRITE(6,5) TITLE 
5 FORMAT(IOAS) 

READ( 50 ,^,END=100 ) M0NTH2 ,STDAY ,YEAR2 
WRITE (6,6) M0NTH2,STDAY,YEAR 
6 FORMAT!///, 

n ■ DATE OF STARTING GPS WEEK (MONTH DAY YEAR) *,13,13,15) 

READS TWO SETS OF PITCH AND ROLL DATA 

READ! 51 ,^,END = 100 ) MONTH!! 1 ) ,DAY1! 1 ) ,YEAR1! 1 ) ,H0UR1! 1 ) ,MIN1! 1 ) , 
$ SECl! 1 ) , PITCH! 1 ) ,ROLL! 1 ) 

TIME1=DFL0AT! DAYK 1 ) )*SECD+DFLOAT( HOURl! 1 ) )*SECH 
$ + DFLOAT!MIN1!1))*60.DO+DBLE!SEC1!1)) 

READ! 51 ,^,END = 100 ) MONTH!! 2 ) ,DAY1! 2 ) ,YEAR1! 2 ) , HOURl! 2 ) ,MIN1! 2 ) , 
$ SECl! 2 ), PITCH! 2 ) ,ROLL! 2 ) 

TIME2=DFL0AT! DAYK 2 ) )*SECD+DFLOAT! HOURl! 2 ) )^SECH 
$ + DFLOAT!MIN1!2))*60.DO+DBLE!SEC1!2)) 

READS TWO SETS OF COURSE 

READ! 52 ,END = 100 ) MONTHS! 1 ) ,DAY3! 1 ) , YEARS! 1 ) , HOURS! 1 ) ,MINS! 1 ) , 

$ SECS! 1), HEAD! 1) 

TIME<^=DFLOAT! DAYS! 1 ) )^SECD+DFLOAT! HOURS! 1 ) )^^SECH 
$ + DFLOAT!MINS( 1) )*60.D0+DBLE!SECS! 1)) 

READ! 52 ,END=100 ) MONTHS! 2 ) ,DAYS( 2 ) , YEARS! 2 ) , HOURS! 2 ) ,MINS! 2 ) , 

$ SECS! 2), HEAD! 2) 

TIME5=DFL0AT! DAYS! 2 ) )^SECD+DFLOAT! HOURS! 2 ) )^SECH 
$ + DFLOAT!MINS!2))^60.D0+DBLE!SECS!2)) 

READS ONE POSITION OF THE SHIP 

20 READ!50,END=100) TT , X ,Y,Z, DLAT , DLOG , DHT , NSAT 

COMPUTES THE DATE 

CALL SUB! TT ,DAY2 ,HOUR2 ,MIN2,SEC2 ) 

DAY2=DAY2+STDAY 

CONVERTS LATITUDE TO DEGREES MINUTES AND SECONDS 

CALL DMS!DLAT,LADEG,LAMIN,LASEG) 

CONVERTS LONGITUDE TO DEGREES MINUTES AND SECONDS 
IN A RANGE 0 TO 180 
IF (DLOG. GT. 180. DO) DLOG=S60 . DO-DLOG 

CALL DMS(DLOG,LODEG,LOMlN,LOSEG) 

WRITE (6,2^) LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG 
2^ FORMAT! ' LAT, LON (DEG MIN SEC, DEG MIN SEG ) ' ,2! 2I^f ,1X,F7.4 ) ,/ ) 

TIMES^DFLOAT! DAY2 )*SECD+DFLOAT! HOUR2 )*SECH 
$ ♦ DFLOAT!MIN2)*60.D0+DBLE!SEC2) 

POSITION IN TIME BETWEEN THE TWO SETS OF PITCH AND ROLL DATA 
AND THE TWO SETS OF COURSE DATA 
25 CONTINUE 

POSITION TO EARLY FOR PITCH AND ROLL DATA 
IF! TIMES. GE. TIMED GO TO 40 

WRITE ! 6 ,S0 ) MONTH2 ,DAY2 ,YEAR2 ,HOUR2 ,MIN2 ,SEC2 
SO FORMAT!* PSN 3 * ,12, IX, 12, IX, 14, IX, 12, IX, 12, IX, F4. 1, 

$ ' REJECTED. NO PITCH AND ROLL DATA') 

GOTO 20 

POSITION TO EARLY FOR PITCH AND ROLL DATA 
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^►0 CONTINUE 

IF(TIME3.GE. TIME‘S) GO TO 60 
HRITE( 6 ,50 ) M0NTH2 ,DAY2 ,YEAR2 ,H0UR2 ,MIN2 ,SEC2 
50 FORMAT ( • PSN 3 M2, IX, 12, IX, IX, 12, IX, 12, IX, F4, 1, 

$ ■ REJECTED. NO COURSE DATA') 

GOTO 20 

POSITION LATER THAN THE SECOND SET OF PITCH AND ROLL DATA 
60 CONTINUE 

IF(TIME3.LE.TIME2) GOTO 70 

MOVES THE PITCH AND ROLL DATA TO THE FIRST SET 

MONTHl(l)=MONTHl(2) 

DAYK 1)=DAY1( 2 ) 

YEAR1(1)=YEAR1(2) 

H0UR1(1)=H0UR1(2) 

MIN1(1)=MIN1(2) 

SEC1(2)=SEC1(2) 

PITCH(1)=PITCH(2) 

ROLL(l)=ROLL( 2) 

TIME1=TIME2 

READS A SECOND SET OF PITCH AND ROLL DATA 

READ( 51,^^,END = 100 ) MONTHK 2 ) ,DAY1( 2 ) ,YEAR1( 2 ) ,H0UR1( 2 ) ,MIN1( 2 ) , 
$ SECK 2 ) ,PITCH( 2 ) ,ROLL( 2 ) 

TIME2=DFLOAT( DAYK 2 ) )i^SECD+DFLOAT( HOURK 2 ) l^SECH 
$ + DFLOAT(MINl( 2 ) )*60.D0+DBLE(SECK 2 ) ) 

GOTO 25 



POSITION LATER THAN THE SECOND SET OF COURSE DATA 
70 CONTINUE 

IF(TIME3.LE. TIMES) GOTO 80 

MOVES THE COURSE DATA TO THE FIRST SET 

MONTH3(l)=MONTH3( 2) 

DAY3(1)=DAY3(2) 

YEAR3(1)=YEAR3(2) 

HOUR3(l)=HOUR3(2) 

MIN3(1)=MIN3(2) 

SEC3(2)=SEC3(2) 

HEAD(1)=HEAD(2) 

TIME4=TIME5 

READS A SECOND SET OF PITCH AND ROLL DATA 

READ( 52 ,END=100 ) MONTH3( 2 ) ,DAY3( 2 ) ,YEAR3( 2 ) ,HOUR3( 2 ) ,MIN3( 2 ) , 
$ SEC3(2),HEAD(2) 

TIME5=DFL0AT( DAY3( 2 ) )*SECD+DFLOAT( HOUR3( 2 ) )^SECH 
$ ♦ DFLOAT(MIN3(2))5^60.D0+DBLE(SEC3(2)) 

GOTO 25 

POSITION HITHIN THE THO SETS IN TIME 

80 RATE=(TIME3-TIME1)/(TIME2-TIME1) 

PIT=PITCH( 1 ) + ( PITCH( 2 )-PITCH( 1 ) )^^RATE 
ROL=ROLL( 1 ) + (ROLL( 2 )-ROLL( 1 ) )^^RATE 
RATE = ( TIME3-TIME^ )/( TIME5-TIME^ ) 

COURSE=HEAD( 1 ) + ( HEAD( 2 )-HEAD( 1 ) )^^RATE 

COMPUTS THE CORRECTIONS AND THE POSITIONS OF THE TRANSDUCER 

CALL PTH(OFFX,OFFY,OFFZ, COURSE, PIT, ROL,DX,DY,DZ) 

COMPUTES THE UTM COORDINATES 
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c 



CALL GPLTTMC LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ,YPOS ,XPOS ) 



XPOS=XPOS+DX 

YPOS=YPOS+DY 

CONVERTS BACK TO GP 

CALL LTTMGP ( YPOS ,XPOS , LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ) 
OUTPUTS THE RESULT 

WRITE ( 7 ) MONTH2 ,DAY2 ,YEAR2 ,HOUR2 ,MIN2 ,SEC2 , 

« LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ,NSAT 

WRITE (6,90) MONTH2 , DAY2 ,YEAR2 , HOUR2 ,MIN2 ,SEC2 , 

« LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ,NSAT 

90 FORMAT ( 2( 1X,I2 ) ,1X,I‘^,2( 1X,I2 ) ,1X,F‘^. 1, 2( 1X,I3 ,1X,I2,1X,F7.4) ,13 ) 

READS ONE MORE POSITION 
GOTO 20 

END OF PROGRAM 

100 CONTINUE 
STOP 
END 

SUBROUTINE PTH( XOFF ,YOFF ,ZOFF , COURSE , PITCH , ROLL ,DX ,DY ,DZ ) 

THIS SUBROUTINE COMPUTES THE CORRECTIONS TO THE 
COORDINATES DUE TO THE OFSET OF ANTENA IN RELATION 
TO THE TRANSDUCER 

PI=ARCOS(-1.0) 

ROLL=ROLL»PI/180. 

PITCH=PITCH»PI/180. 

COURSE=C0URSEitPI/180 . 

CHANGE SIGNAL OF ORIGINAL DATA TO MATCH THE RIGHT HAND 
SYSTEM CONVECTION 

ROLL=-ROLL 

COSPHI= COS(ROLL) 

SINPHI= SIN(ROLL) 

COSK= COS( COURSE) 

SINK= SIN (COURSE) 

COSW= COS( PITCH) 

SINW= SIN( PITCH) 



DX 


= XOFF 


it COSPHI 


it COSK 


* 










$ 


YOFF 


it ( COSW 


it SINK 


+ SINW 


it SINPHI 




COSK 


) + 


$ 


ZOFF 


it ( SINW 


it SINK 


- COSW 


it SINPHI 


* 


COSK 


) 


DY 


» XOFF 


it COSPHI 


it SINK 


it (-1. 


) + 








$ 


YOFF 


it ( COSW 


it COSK 


- SINW 


it SINPHI 


it 


SINK 


) + 


$ 


ZOFF 


it ( SINW 


it COSK 


+ COSW 


it SINPHI 


it 


SINK 


) 



DZ = XOFF it SINPHI - 
$ YOFF ^ ( SINW it COSPHI ) + 
$ ZOFF » ( COSW it COSPHI ) 
RETURN 
END 



C 



SUBROLTTINE GPUTM( LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG , NORTH , E AST ) 
DOUBLE PRECISION A,R,N ,AP ,BP ,CP ,DP ,EP ,S,R1 ,ESQ ,ESQP ,RM,RP ,KO 
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OOU3LE PRECISION R2 ,R5 ,R^,R5,P,P2,PS,P^,P5,P6,A6,B5,SINSEC 
DOUBLE PRECISION PHI ,OLAM,NORTH ,EAST ,B ,PHIMIN ,PI , LON 
REAL LASEG,LOSEG 

THIS SLSROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 125 00 00 W 

CM=-123.000 

PHI=OFLOATI LAOEG )+DFLOAT( LAMIN )/60 . 00+0BLE( LASEG )/3600 . 00 
LON=OFLOAT( LOOEG ) + DFLOAT( LOMIN )/60 . OO+OBLE ( LOSEG )/3600 . DO 
LON =- LON 

DLAM=( LON-CM )*3600 . 000 

A=6378135.0D0 

R=298.2600 

K0=0. 999600 

B = A*(R-1.00)/R 

N = (A-B)/(A+B) 

AP = A*( ( 1.00 -N ) +5. 00/4, 00*(N^H^2-N^«t3 ) +81. 00/64. D0*( N^HH:f-N**S ) ) 

BP = 3. 00/2. 00^^A*((N-N#*2) + 7. 00/8. 00*(N**3-N^)+55. 00/64. D0*N^H^5) 
CP = 15. 00/16. 00*A*(N**2-N**3+3. 00/4. 00*( N**4-N^«f5 ) ) 

OP = 35. 00/48. 00*A*(N**3-N**4+11. 00/16. 00*N**5) 

EP = 315. 00/512. D0^A*(N**4-N**5) 

PHIMIN = PHI*60. 00*2. 908882086660-4 

PI=OARCOS( -1.00 ) 

PHI=PHI/180.00*PI 

S = AP*PHIMIN-BP*OSIN( 2.00*PHI)+CP*DSIN(4.D0*PHI) 

$ - DP*DSIN( 6.00*PHI )+EP*0SIN(8.00*PHI ) 

R1 = KO*S 

SINSEC = ( 1.00/3600.00 )/180.00*PI 
SINSEC=OSIN(SINSEC ) 

ESQ = ( A**2-B**2 )/A**2 
ESQP = ESQ/(1.00-ESQ) 

RM = A*(l. DO-ESQ )/(DSQRT(l.D0-ESQ*DSIN(PHI)**2))**3 
RP = RM*(1.00+ESQP*0COS(PHI)**2) 

R2 = RP*DSIN(PHI)*0COS(PHI)*SINSEC**2/2.00*KO*1.008 
R3 = SINSEC**4*RP*OSIN(PHI)*DCOS(PHI)**3/24.DO*(5.DO-DTAN(PHI)**2 
$ + 9.DO*ESQP*DCOS( PHl)**2+4.00*ESQP*ESQP*0COS( PHI )**4)*K0*1.D16 
R4 = RP*DCOS(PHI )*SINSEC*K0*1.D4 

R5 = SINSEC**3*RP*0COS(PHI)**3/6.00*(1.00-DTAN(PHI)**2 
$ + ESQP*DCOS(PHI )**2)*KO*l.D12 

P = .000100*0LAM 
P2 = P**2 
P3 = P**3 
P4 = P**4 
P5 = P**5 
P6 = P**6 

A6 = P6*SINSEC**6*RP*0SIN(PHI)*0COS(PHI)**5/720.00 
$ * (61. 00-58. D0*0TAN(PHI)**2+DTAN(PHI)**4 

$ + 270 . DO*ESQP*DCOS( PHI )**2-330 . DO*ESQP*DSIN( PHI )**2 )*KO*l . 024 

B5 = P5*SINSEC**5*RP*OCOS( PHI )**5/120. 00*15.00- 18. 00*0TAN( PHI )**2 
$ + 0TAN(PHI)**4+14.00*ESQP*0COS(PHI )**2 

$ - 53.00*ESQP*0SIN(PHI )**2)*KO*1.D20 

N0RTH=R1+R2*P2+R3*P4+A6 
EAST=( R4*P+R5*P3+B5 )+500000 . 00 
RETURN 
END 

SUBROUTINE UTMGP ( NORTH , EAST , LAOEG , LAMIN , LASEG , LOOEG , LOMIN , LOSEG ) 
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DOUBLE PRECISION A,R,N,B, BP ,CP ,DP ,EP ,S ,R1 ,ESQ ,ESQP ,RM ,RP ,K0 
DOUBLE PRECISION R7,R8 ,R9 ,E5,Q,Q2 ,Q5 ,Q4,Q5 ,Q6 ,D6 , AP ,SINSEC 
DOUBLE PRECISION RIO »DLAM, NORTH , EAST ,RPHI ,DPHI ,PHIMIN 
DOUBLE PRECISION EPRIME ,DELTA,D6NUM,CM,PHI ,LON,PI 
REAL LASEG,LOSEG 

THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 125 00 00 H 

CM=-123.0D0 

A=6378135.0D0 

R=298.26D0 

KO=0.9996D0 

PI=DARCOS(-1.DO) 

B = A*(R-1.D0)/R 

N = (A-B)/(A+BJ 

AP = A*( (l.DO-NJ+5. DO/4. D0»(N**2-N**5)+81. DO/64. D0*(N**4-N**5JJ 
BP = 5. DO/2. D0*A^t( (N-N**2J+7. DO/8. D0*(N**3-N**4)+55. DO/64. D0*N**5) 
CP = 15. DO/16. D0*A*(N**2-N*#5+5. DO/4. D0*(N**4-N**5J) 

DP = 35. DO/48. D0*A*(N**5-N**4+11. DO/16. D0*N**5) 

EP = 315.D0/512.D0*A*(N^^*4-N**5J 

FIRST APROXIMATION OF PHI 

PHI =NORTH/30 . 8D0/3600 . DO 

COMPUT TRUE MERIDIONAL DIST AND APROXIMATE PHI 

DO 100 1=1,8 

PHIMIN=PHI*60 . D0*2 . 90888208666D-4 
RPHI=PHI/180.D0*PI 

S = AP*PHIMIN-BP*DSIN(2.D0*RPHI)+CP*DSIN(4.D0*RPHI) 

$ -DP*DSIN( 6 . D0*RPHI ) +EP*DSIN( 8 . DO^^RPHI ) 

R1 = KO*S 
DELTA=N0RTH-R1 
PHI=( DELTA/30.8/3600. )+PHI 
100 CONTINUE 

RPHI=PHI/180.D0*PI 



SINSEC = (l.D0/3600D0)/180.D0#PI 

SINSEC = DSIN( SINSEC J 

ESQ = ( A**2-B**2 J/A**2 

ESQP = ESQ/(1. DO-ESQ 3 

RM = A*(1.D0-ESQ)/(DSQRT(1.D0-ESQ«DSIN(RPHI)**2) 3**3 

RP = RM*(1.DO+ESQP*DCOS(RPHI3**23 

R7=( DTAN( RPHI 3/( 2 . D0*RP**2*SINSEC 3 3*( 1 . DO+ESQP*DCOS( RPHI 3**2 3 
$ *l.D12/KO**2 

R8=( DTAN( RPHI 3/( 24. D0*RP**4*SINSEC 3 3*( 5 . DO+3 . D0*DTAN( RPHI 3**2 
$ +6 . DO*ESQP*DCOS( RPHI 3**2-6 . DO*ESQP*DSIN( RPHI 3**2-3 . D0*ESQP**2 

$ *DCOS( RPHI 3**4-9 . DO*ESQP**2*DCOS( RPHI 3**2*DSIN( RPHI 3**2 3 
$ *1.D24/K0**4 

R9=l . DO/DCOS( RPHI 3/( RP*SINSEC 3*1 . D6/K0 

R10 = l . DO/DCOS( RPHI 3/( 6 . D0*RP**3*SINSEC 3*( 1 . DO+2 . D0*DTAN( RPHI 3**2 
$ +ESQP*DCOS( RPHI 3**2 3*1 . D18/K0**3 

EPRIME=EAST-500000.D0 

Q = .000001D0*EPRIME 

Q2 = Q**2 

Q3 = Q**3 

Q4 = 

Q5 = Q**5 
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Q6 = 

D6NUM=Q6J^DTAN(RPHI ) 

D6 = ( D6NUM )/( 720 . D0*RP**6^^SIN$EC )*( 61 , DO+90 . DO^OTANt RPHI )^^*2 
$ +45.D05^DTAN(RPHI )?^»4+107 . DO»ESQP*OCOS( RPHI 
$ -162.DO^^ESQP5^DSINtRPHI )**2-45. DO^^ESQP^^DTANt RPHI )**2 
$ *OSIN(RPHI )^2 )^^1.D36/K0**6 
E5=( Q5^^1 . 00/0COS( RPHI ) )/( 120 . 00^(RP**5*SINSEC )*( 5. DO+28. DO 
$ 5^DTAN( RPHI )Jt^2 + 2<h. OOJ^DTANt RPHI )^^*4+6 . OO^^ESqP^^DCOSt RPHI 
$ +8 . DO*ESQP^^DSIN( RPHI 2 )^^1 . D50/K0^^*5 

C 

DPHI = ( -R7*Q2+R8*Q4-D6 )/3600.D0 

OLAM=(R9)^Q-R10?^Q3+E5}/3600.00 

PHI=PHI+DPHI 

LON=CM+DLAM 

CALL DMS(PHI,LADEG,LAMIN,LASEG) 

CALL OMS(LON,LOOEG,LOMIN,LOSEG) 

RETURN 

ENO 

C 



C 



SUBROUTINE OMS( DEC ,LOEG,MIN,SEC ) 



C 



DOUBLE PRECISION DEC ,XNUM,XMIN 



XNUM=DABS(DEC) 

LOEG=DINT(XNUM) 

XMIN=(XNUM-DFLOAT( LDEG) )*60.DO 
MIN=DINT(XMIN) 

XNUM=( XMIN-OFLOAT( MIN ) )?^60 . 00 
SEC=SNGUXNUM) 

IF ( LDEG.GE.360) LDEG=LDEG-360 

RETURN 

END 



SUBROUTINE SUB( TT , DAY , HOUR , MIN , SEC ) 

INTEGER DAY, HOUR, MIN 
DOUBLE PRECISION TIME ,TT 
REAL SEC 
TIME=TT 

DAY = DINTtTIME / ( 24 . D0*3600 . DO ) ) 

TIME = TIME - DFLOAT(DAY) * 24, DO ^ 3600. DO 
HOUR = DINK TIME / 3600. DO) 

TIME = TIME - DFLOAT(HOUR) ^ 3600. DO 
MIN = DINT (TIME / 60.00) 

TIME = TIME - DFLOAT(MIN) ^ 60. DO 
SEC = SNGL(TIME) 

RETURN 

END 



//GO.FT06F001 DD 
//GO.FT07F001 DO 
//GO.FT30F001 DD 
//GO.FT50F001 OD 
//GO.FT51F001 DD 
//GO.FT52F001 DD 
// 



SYS0UT=5^ 

DSN=MSS . S0812 . GPS . TRANSDUC . POS , DI$P=SHR 
DSN=MSS . S0812 . CVDATA ,DISP=SHR 
DSN=MSS .S0812 . GPS . ANTEhWA . POS ,DISP=SHR 
DSN=MSS . S0812 . PITROLL ,OISP=SHR 
OSN=MSS . S0812 . SOAS . COURSE , OISP=SHR 
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APPENDIX G 

PROGRAM COMPARE POSITION. SOURCE LISTING 



//JOBCOMPA JOB (0812,9999), ‘E2EQUIEL' ,CLA$S=C 
//*MAIN 0RG=NPGVM1 . 0812P ,LINES=( 99 ) ,CAR0S=( 99 ) 

//^FORMAT PR,0DNAME=GO.FT06F001, 

//*0EST=NPGBACH 
// EXEC FORTVCLG 

//FORT.SYSIN DO * 

C 

C PROGRAM COMPARE POSITION 
C 

C AUTHOR: AUGUSTO EZE(5UIEL 

C GATE: MARCH 27, 1987 

C 

C DESCRIPTION: 

C 

C THIS PROGRAMS TAKES EACH GPS POSITION , FINDS THE 

C CORRESPONDING POSITION USING THE MR FALCON, COMPUTES THE 

C X, Y AND RADIAL DIFFERENCES AND PRINTS THE DATE TIME TAG 

C GEOGRAPHIC POSITIONS AND DIFFERENCES. THIS PROGRAM READS 

C THE PERIODS TO COMPARE THE DATA FROM A SEPARATE FILE. 

C 

C THE PROGRAM RUNS IN MVS 

C 

C I/O SPECIFICATIONS: SEE END OF THIS JOB. 

C 

C THIS PROGRAM IS RESTRICT AS IS TO THE SEAFLOOR BENCHMARK 
C EXPERIMENT 

C 

C THE PROGRAM WILL RUN WITH ANY AMOUNT OF DATA ONLY LIMITED TO DISK 
C SPACE 

C 

C ANY BLANK LINES WILL TERMINATE THE PROGRAM IN ERROR 

C 

c 

c 

DOUBLE PRECISION XPOS( 3 ) ,YPOS( 3 ) ,TIME( 5 ) ,SECD,SECH ,RATE 
DOUBLE PRECISION DIFX,DIFY ,DIST ,XNEW,YNEW 
REAL SEC(5),LASEG(3),L0SEG(3) 

INTEGER MONTHC 5 ) ,DAY( 5 ) ,YEAR( 5 ) ,H0UR( 5 ) ,MIN( 5 ) ,NSAT 
INTEGER LADEG ( 3 ) , LAMIN( 3 ) , LODEG( 3 ) , LOMIN( 3 ) , LINE , LOOP , INDEX 
C 

C INITIALIZATION OF CONSTANTS 
C 

SECH=3600.D0 

SECD=24.D0*SECH 

DIFX=0.0 

DIFY=0.0 

DIST=0.0 

C 

C INPUTS THE NUMBER OF PERIODS TO COMPARE THE DATA 
C 

READ(2S,»,END=100) LOOP 
C 

C LOOP OF THE PROGRAM TO THE NUMBER OF TIMES DESIRED 
C 

DO 1000 INDEX = 1, LOOP 
C 

LINE=60 
REWIND (30) 

REWIND(SO) 

C 

READ( 2S,*,END = 100 ) MONTH( ^ ) ,DAY( ^ ) ,YEAR( ) ,HOUR( 4 ) ,MIN( ^ ) ,SEC( ^ ) 
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« ,MONTH( 5 ) ,DAY( 5 ) ,YEAR( 5 ) ,HOUR( 5 ) ,MIN( S ) ,SEC( 5 ) 

STARTING TIME 

TIME( )=DFLOAT( 0AY( ) )*SECD+OFLOAT( HOUR( ) )»SECH 
S + DFLOATIMINt't) )»6D.DD+0BLE(SEC(<^) ) 

ENDING TIME 

TIME( S )=DFLOAT( DAY! S ) )»SECD+DFLOAT( HOUR( S ) )»SECH 
$ ♦ DFLOAT(MIN(5))»^60.DD+DBLE(SEC(5) ) 

READS TWO MR FALCON POSITIONS 

ID READ(30,END=1DD) MONTHd ),DAY( 1),YEAR( 1),H0UR( 1 ),MIN( 1),SEC( 1), 

$ LADEGI 1 ),LAMIN( 1) >LASE6( 1 ) ,L0DE6( 1 )>L0MIN( 1 )aOSEG( 1) 

TIME( 1 )=DFLOAT( DAY( 1 ) )*SECD+DFLOATl HOUR! 1 ) )»SECH 
$ + DFL0AT(MIN(in»60.D0+DBLE(SEC(D) 

POSITIONS THE FIRST POSITION INSIDE THE PERIOD. THE REST, OF THE 
PROGRAM WILL ADJUST BY IT SELF 

IF( TIME(H.LT.TIME(<i) ) GO TO 10 



READ( 3D ,END=100 ) MONTH! 2 ) ,DAY( 2 ) ,YEAR( 2 ) ,H0UR( 2 ) ,MIN( 2 ) ,SEC( 2 ) , 

$ LADEGI 2 ) ,LAMIN( 2 ) ,LASE6( 2 ) ,L0DE6( 2 ) ,LOMIN( 2 ) ,L0SE6I 2 ) 

TIME! 2 )=DFL0AT( DAY! 2 ) )»SECD+DFLOAT( HOUR! 2 ) )»SECH 
« + DFLOAT(MIN(2))»t60.D0+DBLE(SEC(2n 

COMPUTES THE UTM COORDINATES 

CALL GPUTMI LADEGI 1 ) ,LAMINI 1 ) ,LASEGI 1 ) ,L00EGI 1 ) ,L0MINI 1 ) ,L0SEGI 1 ) , 
« YPOSI D.XPOSI 1 n 

CALL GPUTMI LADEGI 2 ) ,LAMINI 2 ),LASEGI 2 ) ,L0DEGI 2 ) ,LOMINI 2 ),L0SEGI 2 ) , 
$ YP0SI2),XP0SI2)) 

READS GPS POSITION OF THE SHIP 

2D READI SD ,END=1DD ) MONTHI 3 ) ,DAYI 3 ) ,YEARI 3 ) ,H0URI 3 ) .MINI 3 ) ,SECI 3 ) , 

S LADEGI 3 ) , LAMINI 3 ) , LASEGI 3 ) .LODEGI 3 ) .LOMINI 3 ) , LOSEGI 3 ) ,NSAT 

CALL GPLTTMI LADEGI 3 ) , LAMINI 3 ) , LASEGI 3 ) , LODEGI 3 ) , LOMINI 3 ) , LOSEGI 3 ) , 
$ YP0SI3),XP0SI3)) 

TIME! 3 )=DFL0ATI DAYI 3 ) )*SECD+DFLOATI HOURI 3 ) )*SECH 
$ ♦ DFL0ATIMINI3)Mf60.DD+DBLEISECI3)) 

POSITION OUTSIDE OF THE PERIOD 
PROGRAM WILL PICK ANOTHER PERIOD 

IFI TIMEI3).GT.TIMEIS) ) GO TO IDDD 

POSITION IN TIME BETWEEN THE TWO MR FALCON POSITIONS 

25 CONTINUE 

POSITION TO EARLY FOR MR FALCON DATA 
IFITIMEI3 I.GE.TIMEI 1) ) GO TO 40 

* WRITEI6.30) MONTHI3),DAYI3),YEARI3),HOURI3),MINI3),SECI3) 

» 30 FORMAT I ' GPS PSN 3 ' ,12 ,1X,I2 ,1X,I4 ,1X,I2 ,1X,I2 ,1X,F4. 1 , 

» $ ' REJECTED. NO MR FALCON DATA') 

GOTO 2D 

POSITION LATER THAN THE SECOND MR FALCON POSITION 
4D CONTINUE 

IFITIMEI3).LE.TIMEI2)) GOTO 70 
MOVES THE SECOND MR FALCON POSITION TO THE FIRST SET 
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M0NTH(1)=M0NTH(2) 

DAY(1)=DAY(2) 

YEAR(1)=YEAR(2) 

H0UR(1)=H0UR(2) 

MIN(1)=MIN(2) 

SECC2)=SEC(2) 

TIME(1)=TIME(2) 

LADEGI1) = LADEG( 2) 

LAMINC1)=LAMIN(2) 

LASEGI1)=LASEG(2) 

YP0S(1)=YP0S(2) 

LODEG( 1)=LOOEG< 2) 

L0MIN(1)=L0MIN(2) 

L0SEG(1)=L0SEG( 2) 

XP0S(1)=XP0S(2) 

READS A SECOND SET OF MR FALCON DATA 

READ( 30 ,END = 100 ) MONTH! 2 ) ,DAY( 2 ) ,YEAR( 2 ) ,H0UR! 2 ) ,MIN( 2 ) ,SEC( 2 ) , 

$ LADEG! 2 ) , LAMINt 2 ) , LASEGC 2 ) , LODEGI 2 ) , LOMIN! 2 ) , L0SEG( 2 ) 

TIME! 2 )=DFLOAT! DAY! 2 ) )*SECD+DFLOAT! HOUR! 2 ) )*SECH 
$ + 0FLOAT!MIN!2))*60.D0+DBLE!SEC!2)) 

CALL GPUTM! LADEG! 2 ) , LAMIN! 2 ) , LASEG! 2 ) , LODEG! 2 ) , LOMIN! 2 ) , LOSEG! 2 ) , 
$ YPOS!2),XPOS!2)) 

GOTO 25 

GPS POSITION WITHIN THE TWO MR FALCON POSITIONS 

70 RATE = ! TIME! 3 )-TIME! 1 ) )/! TIME! 2 )-TIME! 1 ) ) 

XNEW=XPOS! 1) + ! XPOS! 2 )-XPOS! 1 ) )*RATE 
YNEW=YPOS! 1 ) + ! YPOS! 2 )-YPOS! 1 ) )*RATE 

DIFX=XPOS!3)-XNEW 
DIFY=YPOS!3)-YNEW 
DIST=DSQRT! DIFX?h^2+DIFY^2 ) 

OUTPUTS THE RESULT 



IF !LINE.LT.60) GO TO 85 
WRITE! 6,80) 

80 FORMAT! 1H1,/,17X, 'DIFFERENCES BETWEEN GPS AND MR FALCON POSITIONS' 

$ ,//,29X, 'AVOIDING SVll' , 

$ ,//,4X,'DATE' ,6X, 'TIME' ,09X ,' LATITUDE ' ,7X, ' LONGITUDE ' , 

$ 12X, '! METERS )• ,9X, ,/,2X, 'M' , 2X, ' D ' ,1X, ' YEAR ' ,1X,'HH' , 

$ IX, 'MM' ,1X, 'SS.S' ,2X, 'DD' ,1X, 'MM MX, 'SS.SSSS' ,3X, 'ODD' , IX, 

$ 'MM' ,1X, 'SS.SSSS* ,5X, 'DX' ,5X, ' DY ' ,5X, ' DIST . ' ,1X, 'SAT. * ,/ ) 

LINE=0 

85 LINE=LINE+1 

WRITE! 6 , 90 ) MONTH! 3 ) ,DAY! 3 ) ,YEAR! 3 ) ,HOUR! 3 ) ,MIN! 3 ) ,SEC! 3 ) , 

$ LADEG! 3 ) , LAMIN! 3 ) , LASEG! 3 ) , LODEG! 3 ) , LOMIN! 3 ) , LOSEG! 3 ) , 

$ DIFX,DIFY,DIST,NSAT 

90 FORMAT! 2! IX ,12 ) , 1X,I^, IX ,12 , IX ,12 , IX, F^. 1 ,1X ,13 , IX , 12 ,1X, F7.^f , 

$ • N* ,1X,I3,1X,I2,1X,F7.^, ' W' ,1X,F6.2,1X,F6.2,1X,F7.2,1X,I2) 

READS ONE MORE POSITION 
GOTO 20 

END OF PROGRAM 

1000 CONTINUE 
100 CONTINUE 
STOP 
END 



C 



SUBROUTINE GPUTM! LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG , NORTH , EAST ) 
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DOUBLE PRECISION A ,R ,N, AP ,BP ,CP , DP ,EP ,S ,R1 ,ESQ ,ESQP ,RM ,RP ,KO 
DOUBLE PRECISION R2,R3 ,R4,R5,P,P2,P3 ,P4,P5,P6 ,A6 ,B5,SINSEC 
DOUBLE PRECISION PHI ,DLAM,NORTH ,EAST ,B ,PHIMIN,PI , LON 
INTEGER LADEG,LAMIN,LOOEG,LOMIN 
REAL LASEG,LOSEG 

THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 123 00 00 W 

CM=-123.0D0 

PHI=DFLOAT( LADEG)+DFLOAT( LAMIN )/60 . DO+DBLE ( LASEG 1/3600.00 
LON=DFLOAT( LODEG )+DFLOAT( LOMIN 1/60 . DO+DBLE ( LOSEG 1/3600. DO 
LON=-LON 

DLAM=( LON-CM 1^^3600 . 000 
C 

A=6378135.000 

R=298.26D0 

C 

K0=0. 999600 
C 

B = A^^( R-l.OOl/R 
C 

N = (A-B1/(A+B1 

AP = A*( ( 1.00-N1+5.D0/4.D0*(N5^*2-N^^5^31+81.D0/6^.00*(N**4-N**51 1 
BP = 3 . DO/2 . ( N-N*^( 2 1 +7 . DO/8 . D0*( N^^^^3-N^H^ 1 +55 . 00/6^ . 00*N^^*5 1 

CP = 15. 00/16. D0*A^^(N*i(2-N*^(3+3. DO/4. 00^^( N?«^4-N**5 1 1 
OP = 35. DO/48. 00^^A*(N^^*3-N**4 + 11. DO/16. DO^^NiH^Bl 
EP = 315.00/512.D0^(A*(N^(*4-N*i^51 
PHIMIN = PHI*60. 00*2.908882086660-4 
C 

PI=OARCOS( -1.001 
C 

PHI=PHI/180.D0*PI 

C 

S = AP*PHIMIN-BP*OSIN( 2.00*PHI 1+CP*DSIN(4.D0*PHI 1 
$ - 0P*DSIN(6.D0*PHI1+EP*0SIN(8.00*PHI1 

R1 = KO*S 
C 

SINSEC = ( 1.00/3600. DO 1/180. DO*PI 
SINSEC=DSIN(SINSEC1 
C 

ESQ = ( A**2-B**21/A**2 
ESQP = ESQ/( l.OO-ESQl 

RM = A*(1.D0-ESQ1/(DSQRT(1.D0-ESQ*0SIN(PHI1**211**3 
C 

RP = RM*( 1.00+ESQP*DCOS( PHI 1**21 

R2 = RP*DSIN(PHIl*0COS(PHIl*SINSEC**2/2.00*KO*1.0D8 
R3 = SINSEC**4*RP*OSIN(PHI1*DCOS(PHI1**3/24.DO*(5.DO-DTAN(PHI1**2 
$ + 9.DO*ESQP*OCOS( PHI l**2+4.00*ESQP*ESQP*0COS( PHI 1**4 1*K0*1.D16 
R4 = RP*DCOS( PHI 1*SINSEC*KO*1.04 

R5 = SINSEC**3*RP*OCOS(PHI1**3/6.DO*(1.DO-DTAN(PHI1**2 
$ + ESQP*OCOS( PHI 1**2 1*K0*1. 012 

C 

P = .000100*0LAM 
P2 s p^Z 
P3 = P**3 
P4 = P**4 
P5 = P**5 
P6 = P**6 
C 

A6 = P6*SINSEC**6*RP*0SIN( PHI 1*DC0S( PHI 1**5/720.00 
$ * (61.00-58.00*DTAN(PHI1**2+0TAN(PHI1**4 

$ + 270 . 00*ESQP*0COS ( PHI 1**2-330 . DO*ESQP*DSIN( PHI 1**2 1*K0*1 . D24 

B5 = P5*SINSEC**5*RP*OCOS( PHI 1**5/120 . D0*( 5. DO-18 . 00*DTAN( PHI 1**2 
$ + DTAN( PHI 1**4+14.DO*ESQP*OCOS( PHI 1**2 

$ - 58.00*ESQP*DSIN( PHI 1**2 )*KO*1.D20 

N0RTH=R1+R2*P2+R3*P4+A6 
E AST= ( R4*P+R5*P3+B5 1+500000 . DO 
RETURN 
END 
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END 



//GO. FT06F001 DO SYSOUT=J* 

//G0.FT25F001 DO 0SN=MSS.S0812. PERIODS. COMPARE, DISP=SHR 
//GO.FT3OF0O1 00 OSN=MSS .S0812 . FALCON. TRANSDUC . POS,OISP=SHR 
//GO.FTSOFOOl DO DSN=MSS.S0812. GPS. TRANSDUC. POS,01SP=SHR 
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APPENDIX H 

PROGRAM COMPARE PLOT. SOURCE LISTING 



//E2EQUIEL JOB ( 0812 , 9999 E2EQUIEL SCLASS=C 

// EXEC FRTVCLGP 

//FORT.SYSIN DD * 

C 

C PROGRAM COMPARE PLOT 

C 

C RUNS IN FORTRAN VS 

C 

C AUTHOR; AUGUSTO E2EQUIEL 

C 

C DATE : 31 MARCH 1987 

C 

C THIS PROGRAM MAKES A PLOT OF THE POSITIONS OF THE SHIP 
C USING GPS AND MR FALCON POSITIONS OF THE TRANSDUCER 
C 

c 

DOUBLE PRECISION XPOS,YPOS,XLEFT,YLEFT,TIME( 4 ) ,SECD,SECH 
REAL SEC( 3 ) ,XPLT , YPLT ,BLXH ,BLYH , SCALE , VALUE , LASEG, LOSEG 
REAL YORGjXORG 

INTEGER M0NTH( 3 ) ,YEAR( 3 ) ,0AY( 3 ) ,HOUR( 3 ) ,MIN( 3 ) 

INTEGER IPEN , LADEG , LAMIN , LODEG , LOMIN 
INTEGER 10, LOOP, INDEX 

DIMENSIONS OF SHEET, LEFT CORNER AND SCALES 

BLXH=1S. 

BLYH=20.5 
SCALE =l./5000. 

XORG=0.0 
Y0RG=-25.0 

INITIALIZATION OF CONSTANTS 

TIME(4)=0.D0 
SECH=3600.D0 
SECD=24.D0*SECH 

PLOTTER INITIALIZATION 

CALL PLOTS (0,0,0) 

READS THE AMOUNT OF PERIODS TO COMPARE 
READ(30,») LOOP 
DO 100 INDEX=l,LOOP 

READS THE PERIODS AND INFERIOR LEFT CORNER OF EACH AREA 

READ( 30 ,* ) MONTH( 1 ) ,DAY( 1 ) ,YEAR( 1 ) ,HOUR( 1 ) ,MIN( 1 ) ,SEC( 1 ) 

$ ,MONTH(2),OAY(2),YEAR(2),HOUR(2),MIN(2),SEC(2), XLEFT,YLEFT 

STARTING TIME 

TIME( 1 )=DFL0AT( DAY( 1 ) )*SECD+DFLOAT( HOUR( 1 ) )*SECH 
$ + 0FLOAT(MIN(l))*60.D0+DBLE(SEC(l)) 

ENDING TIME 

TIME( 2 ) = DFLOAT( DAY( 2 ) )*SECD+DFLOAT( HOUR( 2 ) )*SECH 
$ 4. 0FLOAT(MIN(2))*6O.D0+DBLE(SEC(2)) 
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POSITIONS THE PEN IN THE ORIGIN 



YORG=-YORG 

XORG=25.0*( l-( -l)**INDEX)/2 
IF(INDEX.EQ.2) XORG=0 
CALL PLOT (X0RG,Y0RG,-3) 

GRID 



CALL PLOT (0.1, 0.1, 3) 

CALL PLOT (0.1,BLYH-0.1,2) 

CALL PLOT (BLXH-0.1,BLYH-0.1,2) 
CALL PLOT (BLXH-0.1,0.1,2) 

CALL PLOT (0.1, 0.1, 2) 



SCALE 



CALL METER( SCALE) 
TITLE 



CALL SYMBOL(08.1,19.2,0.20,29HSEAFLOOR BENCHMARK EXPERIMENT ,0 ., 29 ) 
CALL SYMBOL(08.1,18.8,0.20,29H PHASE II ,0.,29) 

CALL SYMBOL(08.1,18.4,0.20,29H 16 AUGUST 1986 ,0.,29) 

CALL SYMBOL(08.1,18.0,0.20,29HAVOIDING SVll. NO CONSTRAIN ,0.,29) 
CALL SYMBOL( 01.0,03.375,0.15,3,0. ,-l) 

CALL SYMBOL(1.20,03.3,0.15,29HMR FALCON TRANDUCER POSITIONS,0 . ,29 ) 
CALL SYMB0L( 01.0,03.075,0. 15,11,0. ,-l) 

CALL SYMBOL( 1.20,03.0,0. 15, 29HGPS TRANDUCER POSITIONS ,0.,29) 



SETS SYMBOL FOR MINI RANGER POSITIONS 



ISN=3 



REWINDS THE FILES 

REWIND(31) 

REWIND( 32) 

LOOP FOR THE TYPES OF POSITION 



DO 80 10=31,32 
SETS PEN UP 



IPEN=3 

LOOP WITHIN ONE TYPE OF POSITION 
10 CONTINUE 
READ THE DATA 

READ( IO,END=50 )MONTH( 3 ) ,DAY( 3 ) ,YEAR( 3 ) ,HOUR( 3 ) , 

$ MIN( 3 ) ,SEC ( 3 ) , LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG 

COMPUTES THE UTM COORDINATES 

CALL GPUTM ( LADEG , LAMIN , LASEG , LODEG , LOMIN , LOSEG ,YPOS ,XPOS ) 
COMPUTES THE TIME IN SECS 



TIME ( 3 )=DFLOAT( DAY( 3 ) )^^SECD+DFLOAT( HOUR( 3 ) )i^SECH 
$ + DFLOAT(MIN(3))*60.D0+0BLE(SEC(3)) 
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TESTS IF WITHIN THE PERIOD 

IF(TIME(3).GT.TIME(2)) GOTO 50 
IF(TIME(3).LT.TIME(in GOTO 10 

COMPUTES THE PLOTTER COORDINATES 

XPLT = ( XPOS-XLEFT )*SCALE^^100 . 

YPLT = ( YPOS-YLEFT )*SCALE*100. 

TESTE IF INSIDE AREA 

IFIXPLT.LT.O. ) GOTO 30 
IFIYPLT.LT.O. ) GOTO 30 
IF(XPLT.GT.BLXH) GOTO 30 
IF( YPLT.GT.BLYH) GOTO 30 

TESTS IF POSITIONS ARE AWAY MORE THEN 10 SECS IN TIME 

IF(TIME(3)-TIME(4).GT.10.D0 ) IPEN=3 
TIME(4)=TIME(3) 

PLOTS THE POSITION 

CALL PL0T(XPLT,YPLT,IPEN) 

IPEN=2 

PLOTS THE TIME EVERY MINLTTE 

IF(IFIX(SEC(3)).GT.O) GO TO 10 
VALUE =FLOAT ( HOUR( 3 ) )^^100 + FLOAT( MIN( 3 ) ) 

CALL NUMBER(XPLT+0. 15, YPLT+0. 15, 0.15, VALUE ,0. ,-l ) 
CALL SYMB0L(XPLT,YPLT,0.1S,ISN,0. ,-l) 

GOTO 10 

PEN UP WHILE THE POSITIONS ARE OUT OF THE SHEET 

30 IPEN=3 

GOTO 10 

END OF PLOT 

50 CONTINUE 

ISN=11 
IPEN=3 
80 CONTINUE 
100 CONTINUE 

CALL PLOT! 0. ,0. ,+999) 

STOP 

END 

SUBROUTINE METER( SCALE) 

X0=2.0 

Y0=1.0 

CALL PL0T(X0,Y0,-3) 

DO 10 J=1,12 
XP= FLOAT(J)-1.0 
CALL PLOTC XP,0.0,3) 

CALL PLOT(XP,0.2S,2) 

10 CONTINUE 

DO 20 J=l,9 
XP=FL0AT( J)*0.1 
CALL PL0T(XP,0.0,3) 

CALL PL0T(XP,0.2,2) 
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20 CONTINUE 
C 

CALL PLOT(0. ,0.0,3) 

CALL PLOT( 11. ,0. ,2) 

CALL PLOKO. ,0.2,3) 

CALL PLOT(ll.,.2,2) 

C 

VALUE=1. /SCALE/100. 

CALL NUMBER( -0. 25 ,0 . 27,0 . 15, VALUE ,0 . 0 ,-l 1 
CALL NUMBER! 0. 95,0.27,0.15,0.0,0.0,-1 ) 
VALUE=5. O/SCALE/100 . 

CALL NUMBER! 5. 775, 0.27, 0.15, VALUE, 0. 0,-1) 
VALUE=10 . O/SCALE/100 . 

CALL NUMBER! 10 . 7 ,0 . 27 ,0 . 15 ,VALUE ,0 . 0 , -1 ) 
CALL SYMBOL(5.5,-0.25, .15,6HMETERS,0.0,6) 
C 

CALL SYMBOL! 10.0,12.0,1.5,62,0.0,-1) 

CALL SYMBOL! 10.1,12.4,0.2,85,0.0,-1) 

XO=-XO 

YO=-YO 

CALL PL0T!X0,Y0,-3) 

C 

RETURN 

END 



SUBROUTINE GPUTM ! LADE6 , LAMIN , LASEG , LODEG , LOMIN , LOSEG , NORTH , E AST ) 

DOUBLE PRECISION A ,R,N,AP ,BP ,CP ,DP ,EP ,S,R1 ,ESQ,ESQP ,RM,RP ,KO 
DOUBLE PRECISION R2 ,R3 ,R4, R5,P ,P2 ,P3 ,P4 ,P5 ,P6 ,A6 ,B5 ,SINSEC 
DOUBLE PRECISION PHI ,DLAM, NORTH, EAST, B,PHIMIN, PI , LON 
REAL LASEG, LOSEG 

THIS SUBROUTINE COMPUTS THE UTM COORDINATES OF GP IN WGS 72 
IN ZONE 10 CENTRAL MERIDIAN 123 00 00 H 



C 



C 

C 

C 



C 

c 

c 



c 



c 



CM=-123.0D0 

PHI=DFLOAT! LADEG )+DFLOAT! LAMIN )/60 . DO+DBLE! LASEG )/3600 . DO 
LON=DFLOAT! LODEG )+DFLOAT! LOMIN )/60 . DO+DBLE! LOSEG )/3600 . DO 
LON=-LON 

DLAM=! LON-CM )^^3600 . ODO 

A=6378135.0D0 

R=298.26D0 

KO=0.9996D0 

B = A*!R-1.D0)/R 

N = !A-B)/!A+B) 

AP = A^!!l.D0-N)+5. DO/4. D0*!N^2-NJH^3)+81. DO/64. D0*!N**4-N^HfS)) 

BP = 3. DO/2. D0*A*! !N-N^«^2)+7. DO/8. D0^t!N^t*3-NiH^4)+55. DO/64. D0*N^H^5) 
CP = 15.D0/16.D0*A^t!N**2-N^^*3+3.D0/4.D0^^!N**4-N**S) ) 

DP = 35. DO/48. D0^^A^^!N^3-N**4+11. DO/16. D0^^N**5) 

EP = 315.D0/512.D0*Ai^!Ni^*4-N**5) 

PHIMIN = PHI*60.D0*2.90888208666D-4 

PI=DARCOS(-1.DO) 

PHI=PHI/180.D0i*PI 

S = AP*PHlMIN-BP^^DSIN!2.D0^tPHI)+CP^^DSIN!4.D0^^PHI) 

$ - DPJ^DSIN!6.D0*PHI)+EP^^DSIN!8.D0*PHI) 

R1 = KQ^S 

SINSEC = ! l.D0/3600.D0)/180.D0i^PI 
SINSEC=DSIN!SINSEC) 



ESQ = ! A**2-B»*2 )/A**2 
ESQP = ESQ/! 1. DO-ESQ) 



112 



RM = A*( 1. DO-ESQ )/(DSQRT(l. DO-ESQ*DSIN( PHDJh^Z) 

C 

RP = RM*(l.DO+ESQP*DCOS(PHI)**2) 

R2 = RP*DSIN(PHI)i^DCOS(PHI)»SINSEC^^*2/2.DO*KO*1.0D8 
R3 = SINSEC»*<^*RP#OSIN( PHI )*DCOS( PHI )J^*3/2^ . DO*( 5. DO-DTAN( PHI )**2 
$ + 9 . DO*eSQPi^DCOS( PHI )**2+<^ . DO^^ESQP^^ESQP*DCOS( PHI )*KO*l . D16 

s RP*DCOS(PHI)»SINSECi^KO*l.D<h 
RS = SINSEC**3*RP*DCOS(PHI)**3/6.00if(l,00-DTAN(PHI)**2 
$ + ESQP*DCOS(PHI )**2)*KO*l.D12 

C 

P = .0001D0*DUM 
P2 = PJ^*2 
P3 = P**3 

P5 = P**5 
P6 = P^^*6 
C 

A6 = P6i^SINSECJH^6*RP*DSIN(PHI)*DCOS(PHI)^^*5/720,D0 
$ * ( 61 . DO-58 . DO*DTAN( PHI )J^i^2+DT AN( PHI 

$ + 270 . DO^^ESQP*DCOS( PHI )**2-330 , DO*ESQP#DSIN( PHI )**2 )»K0*1 . D2<h 

85 = P5»SINSEC**5*RP*0COS(PHI)i^*5/120.D0i^(5.D0-18.D0»DTAN(PHI)^2 
$ + DTAN(PHI)**4 + 14,D0*ESQP^^DC0S(PHI)**2 

$ - 58.D0*ESQP*DSIN(PHI)**2)*K0*1.D20 

NORTH=Rl+R2i^P2+R3*P4+A6 
EAST=( R^»P+R5*P3+B5) +500000. DO 
RETURN 
END 



//GO.PLOTPARM DD * 

iPLOT XMIN=0. ,XMAX=300. ,YMIN=0. ,YMAX=51. , SCALE =1. ,UNITS=2. 540 SEND 
//G0.FT06F001 DD SYSOUT=* 

//G0.FT30F001 DD DSN=MSS . S0812 . PERIODS . COMPARE ,DISP=SHR 
//GO.FT31F001 DD DSN=MSS . S0812 . FALCON . TRANSDUC . POS ,DISP=SHR 
//GO.FT32F001 DD DSN=MSS . S0812 . GPS . TRANSDUC . POS ,DISP=SHR 
// 
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APPENDIX I 

CONTROL FILE FOR KALMN2 PROGRAM. AS DEFINED BY 

BRO\VN(1986) 

585590. 700000. 16 10 4 1 1 000 

0 -2740.506 -4341.136 3772.0820 .8738E1 0.0 -3.0 0.0 

6378.135 298.26 

1.0 1.0 1.0 l.E-12 2.1 l.E-14 l.E-50 l.E-50 l.E-50 3.0 4.0 5.0 
1010000100 
0.000005 0.0 

110000101 
980. 15. 75. 10. 

0.75E-4 .75E-4 .75E-4 .lE-4 1.0 1.0 1.0 2.0 .lE-03 .lE-04 

1.0 1.0 1.0 1.0 1.0 1.0 l.E-19 l.E-19 l.E-19 

-.031 
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APPENDIX J 

CONTROL FILE FOR CVFICA PROGRAM. 



The first two lines may have a title up to 80 columns. The third line has the 
month day and year of the starting day of the GPS week. The fourth line has the 
number of satellites to be used, and other two constants are kept at zero. The last line 
has the estimate position in WGS 72 coordinates. 

The next is the one used for this thesis; 

SEAFLOOR BENCH MARK PHASE II 
R V POINT SUR 16 AUGUST 1986 

8 10 1986 
4 0. 0. 

-2739.58 -4340.77 3773.14 
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APPENDIX K 

CONTROL FILE FOR COMPARE POSITION AND PLOT PROGRAMS 



The first lines has the number of periods to compare the data. This line must be 
followed by the same number of lines as periods are to compare. Each line the date 
time tags from the start and of the periods has the left inferior corner of the area 
corresponding to each period. 

The next is the one used for this thesis: 



11 



8 


16 


1986 


19 


14 


O 

o 


8 


16 


1986 


OO 


16 


1986 


19 


20 


0.0 


8 


16 


1986 


OO 


16 


1986 


19 


25 


0.0 


8 


16 


1986 


OO 


16 


1986 


19 


30 


0.0 


8 


16 


1986 


CO 


16 


1986 


19 


35 


0.0 


8 


16 


1986 


OO 


16 


1986 


19 


47 


0.0 


8 


16 


1986 


OO 


16 


1986 


19 


57 


0.0 


8 


16 


1986 


OO 


16 


1986 


20 


07 


0.0 


8 


16 


1986 


8 


16 


1986 


20 


13 


0.0 


8 


16 


1986 


CO 


16 


1986 


20 


20 


0.0 


8 


16 


1986 


CO 


16 


1986 


20 


30 


0.0 


8 


16 


1986 



19 20 3.0 567320.0 4038400.0 
19 25 3.0 566920.0 4039000.0 
19 30 3.0 566670.0 4039600.0 
19 35 3.0 566270.0 4040200.0 
19 40 3.0 566020.0 4040800.0 

19 53 3.0 564920.0 4040950.0 

20 05 3.0 564970.0 4040200.0 
20 13 3.0 565570.0 4039500.0 
20 20 3.0 566020.0 4039000.0 
20 25 3.0 566420.0 4038500.0 
20 37 3.0 565985.0 4038400.0 
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